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The magnetic susceptibihty x*(i) a-nd specific heat C{t) versus temperature t of the spin S — 1/2 
antiferromagnetic (AF) alternating-exchange (Ji and J2) Heisenberg chain are studied for the entire 
range < a < 1 of the alternation parameter a = J2/J1 {Ji, J2 > 0, J2 < Ji, t = k^T/J\, 
X* ~ xJ'i-l^g'^ A)- For the uniform chain (a — 1), the high-accuracy S'nd C{t) Bethe ansatz 

data of Kliimper and Johnston (unpublished) are shown to agree very well at low t with the respective 
exact theoretical low-t logarithmic correction predictions of Lukyanov, Nucl. Phys. B 522, 533 
(1998). Accurate (~ 10"*^) independent empirical fits to the respective data are obtained over t 
ranges spanning 25 orders of magnitude, 5 x 10~^^ < i < 5, which contain extrapolations to the 
respective exact t = Q limits. The infinite temperature entropy calculated using our C{t) fit function 
is within 8 parts in 10* of the exact value In 2. Quantum Monte Carlo (QMC) simulations and 
transfer- matrix density-matrix renormalization group (TMRG) calculations of x*{a, t) are presented 
for 0.002 <t<lQ and 0.05 < a < 1, and an accurate (2 x 10~*) two-dimensional (a,t) fit to the 
combined data is obtained for 0.01 < t < 10 and < a < 1. From the \ow-t TMRG data, the spin 
gap A(a) is extracted for 0.8 < a < 0.995 and compared with previous results, and a fit function 
is formulated for < a < 1 by combining these data with literature data. We infer from our data 
that the asymptotic critical regime near the uniform chain limit is only entered for a > 0.99. We 
examine in detail the theoretical predictions of Bulaevskii, Sov. Phys. Solid State 11, 921 (1969), 
for x*{(^y^) a-nd compare them with our results. To illustrate the application and utility of our 
theoretical results, we model our experimental x{T) and specific heat Cp{T) data for NaV205 single 
crystals in detail. The x{T) data above the spin dimerization temperature Tc « 34K are not in 
quantitative agreement with the prediction for the S = 1/2 uniform Heisenberg chain, but can 
be explained if there is a moderate ferromagnetic interchain coupling and/or if J changes with T. 
Fitting the x{T) data using our x'i'^y't) fit function, we obtain the sample-dependent spin gap and 
range A(r = 0)/kB = 103(2) K, alternation parameter 5{0) = (1 - a)/(l -\- a) = 0.034(6) and 
average exchange constant J(0)/A:b = 640(80) K. The S{T) and A(r) are derived from the data. A 
spin pseudogap with magnitude ~ 0.4 A(0) is consistently found just above Tc, which decreases with 
increasing temperature. From our Cp{T) measurements on two crystals, we infer that the magnetic 
specific heat at low temperatures T < 15 K is too small to be resolved experimentally, and that the 
spin entropy at Tc is too small to account for the entropy of the transition. A quantitative analysis 
indicates that at Tc , at least 77 % of the entropy change due to the transition at Tc and associated 
order parameter fluctuations arise from the lattice and/or charge degrees of freedom and less than 
23 % from the spin degrees of freedom. 

PACS numbers: 75.40. Cx, 75.20.Ck, 75.10.Jm, 75.50.Ee 
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I. INTRODUCTION 

An antiferromagnetic alternating-exchange Heisenberg 
chain is one in which nearest-neighbor spins in the chain 
interact via a Heisenberg interaction, but with two anti- 
ferromagnetic ( AF) exchange constants J2 < Ji , Ji , J2 > 
which alternate from bond to bond along the chain; 
the alternation parameter is a = JijJx- Here we will 
be concerned with the magnetic susceptibility x ^-nd 
specific heat C versus temperature T of alternating- 
exchange chains consisting of spins S — 1/2. The uni- 
form AF Heisenberg chain is one limit of the alternat- 
ing chain in which the two exchange constants are equal 
(a = 1, Ji = J2 = J). At the other limit is the iso- 
lated dimer in which one of the exchange constants is 
zero [a = 0). The present work is a combined theo- 
retical and experimental study of x{T) and C[T) of the 
5=1/2 AF alternating-exchange chain over the entire 
range < a < 1 of the alternation parameter, with the 
emphasis on the regime a < 1 at and close to the uniform 
chain limit. This latter regime is relevant for compounds 
showing second order spin dimerization transitions with 
decreasing T. The present work was originally motivated 
by our desire to accurately extract the temperature de- 
pendent energy gap A(r) for magnetic excitations, the 
"spin gap" , from experimental x(T) data for the 5* = 1/2 
chain/ two-leg ladder compound NaV2 05 below its spin 
dimerization temperature Tc ~ 34 K. We found that ex- 
isting theory for the alternating-exchange chain was in- 
sufficient to accomplish this goal. In the present work we 
critically examine the predictions of previous theory, per- 
form the required additional theoretical calculations, and 
then apply the results to extract A(T) at T ^ Tc from our 
x(r) data for NaV2 05 single crystals. We have extended 
the original goal so that we also include theoretical and 
experimental studies of C(T) and how this quantity re- 
lates to x(T). In the remainder of this introduction we 
briefly review the prior theoretical results pertaining to 
x(T) and C(T) of the uniform and alternating-exchange 
chain to place our work in the proper context. We then 
review the experimental and theoretical background on 
NaV2 05 and describe the plan for the rest of the paper. 

A. Theory 

The x{T) and C{T) of both hmits of the S* = 1/2 
AF alternating-exchange Heisenberg chain are known ex- 
actly. For the dimer, the x(T) is given by the exact 
Eq. below and the exact C{T) is also easily calcu- 
lated. The x{T) and C(T) of the uniform chain for T > 
QAJ/k-Q (fee is Boltzmann's constant) were estimated 
from calculations foE chains with < 11 spins by Bonner 
and Fisher in 1964p they extended their results by ex- 
trapolating to r = 0, and in the case of x(T) to the exact 
T = value.EI The exact solution for x(^) of the uniform 
chain was obtained using the Bethe ansatz in 1994 by 



Eggert, Affleck and Takahashi, and compared with their 
low-T results from conformal field theory.B They found, 
remarkably, that x(T ^ 0) has infinite slope. Their nu- 
merical x{T) values are up to ~ 10% larger than the 
Bonner-Fisher extrapolation for T < 0.25J/A;b (for a 
comparison of the two predictions, see Fig. 8.1 in Ref.|). 
Their conformal field theory calculations showed that the 
leading order correction to the zero temperature limit 
is of the form x{T) = x(0){l + l/[21n(To/T)]}, where 
the value of the scaling temperature Tq is not predicted 
by the field theory. Such log terms are called "logarith- 
mic corrections" in the literature. One of us recently 
presented numerical Bethe ansatz calculations of x{T) 
with a higher absolute accuracy for x(^) estimated to 
be 1 X 10"'' J3 and showed that the data are consistent 
with the above field theory prediction, with an additional 
higher order logarithmic correction, over the temperature 
range 5 x lO'^^ < UbT / J < lO'^. Corresponding C{T) 
calculations were also carried out, and logacithmic cor- 
rections were studied for this quantity as welLtl Lukyanov 
has recently presented an exact theory foij-|\^(T) and C{T) 
at low T, including the exact value of TqB In the present 
work, we compare the very recent-Jiumerical Bethe ansatz 
results of Kliimper and JohnstonB with the predictions of 
Lukyanov's theory and find agreement for xC^) to high 
accuracy (< 1 x 10"^) over a temperature range spanning 
18 orders of magnitude, 5 x lO^^s < ksT/J < 5 x lO"'^; 
the agreement in the lower part of this temperatures 
range is much better, 0(10-'^). For C(T), the logarith- 
mic correction in Lukyanov's theory is insufficient to de- 
scribe the Bethe ansatz data sufficiently accurately even 
at very low temperatures, so we derive the next two log- 
arithmic corrections from the Bethe ansatz C(T) data. 
For various applications, it would be desirable to have 
fits to the x{T) and C(T) Bethe ansatz data which ex- 
tend to higher temperatures. We describe the formula- 
tion and implementation of fit functions, incorporating 
the influence of the logarithmic correction terms, which 
yield extremely precise fits to the data for both quan- 
tities over the entire 25 decades in temperature of the 
calculations, 5 x 10"^^ < ksT/J < 5. 

The x(^) in the intermediate regime < a < 1 
has been investigated analytically in the Hartre&-£'ock 
approximationll and using numerical techniques .u'EI Of 
particular interest here is the regime a 1, close to 
the uniform limit, which is the regime relevant to mate- 
rials exhibiting a dimerization transition with decreasing 
T such as occurs in materials exhibiting a spin-Peierls 
transition. There are no accurate theoretical predictions 
available for x{T) of the alternating-exchange Heisenberg 
chain in this regime, which is the property usually used to 
initially characterize the occurrence of such a transition 
experimentally. To address this deficiency and to also 
cover a more extended a range, we carried out extensive 
quantum Monte Carlo (QMC) simulations and transfer- 
matrix density-matrix renormalization group (TMRG) 
calculationsBlB of xiT) for 0.05 < a < 1 over the tem- 
perature range 0.002 < ksT/Ji < 10. 
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An interesting issue is how the spin gap A evolves 
with alternation parameter a as the uniform limit is ap- 
proached, a ^ 1. Because the uniform chain is a gapless 
quantum-critical system, the introduction of alternating 
exchange along the chain has been theoretically predicted 
to yield a nonanalytic A (a) behavior for a ^ 1. We de- 
rive A(a) by fitting our low-i TMRG xi^) data by an 
expression which we formulated. The A(a) results are 
compared with those of previous numerical calculations 
and with the theoretical prediction. We infer from our 
data that the asymptotic critical regime is only entered 
for a > 0.99. 

In order to be optimally useful for accurately modeling 
experimental x(^) data for alternating-exchange chain 
compounds, our QMC and TMRG x(q!, T) results must 
first be accurately fitted by a continuous function of both 
a and T. We will introduce a general fit function which 
eventually proves capable of fitting these combined data 
for the alternating-exchange Heisenberg chain very accu- 
rately. We first fit the x(^) of the uniform chain and 
isolated dimer using this function and then use the ob- 
tained fitting parameters as end-point parameters in the 
fit to our combined QMC and TMRG data for interme- 
diate values of a. The final fit function is a single two- 
dimensional function of a and T for < a < 1 which can 
be used to extract the (possibly temperature-dependent) 
alternation parameter, exchange constants and spin gap 
from experimental x(T) data for compounds for which 
the S = 1/2 AF alternating-exchange Heisenberg chain 
Hamiltonian is appropriate. Our fit function will also be 
useful as a reference for x(r) calculated from other re- 
lated S = 1/2 Hamiltonians such as that incorporating 
the spin-phonon interaction for spin-Pcierls systems. 

B. NaVaOs 

Vanadium oxides show a remarkable variety of elec- 
tronic behaviors. For example, the metallic fee normal- 
spinel structure compound LiV204 shows local moment- 
like behaviors above ~ 50 K, C|Cassing over to heavy 
fermion behaviors below ^ lOK.EHl On the other hand, 
the compound CaV205 has a two-leg trellis-ladder- 
layer structureE^ in which all of the V atoms are crystallo- 
graphically equivalent and is a Mott-Hubbard insulator. 
The x{T) shows a spin-gap A//cb ~ 660 K arising from 
strong coupling of the Y S — 1/2 spins across a rung.t^ 
Modeling of x(T) by QMC simulations confirmed that 
this compound consists magnetically of V2 dimers, with 
an intradimer AF exchange constant JfkB ~ 665 K and 
with very weak interdimer interactions. 

The compound NaV205 can also be formed. The 
crystal structure was initially reported in 1975 to 
consist of two-leg ladders as in CaV205, but in a 
non-centrosymmetric (acentric) structure (space group 
P2i mn) in which charge segregation occurs such that one 
leg of each ladder consists of V^"* and the other of crystal- 



lographically inequivalent V"*"^ ions.ll3 However, recently 
five different crystal structure investigations showed that 
the structure is actually centrosymmetric (space group 
Pmmn), with all V atoHis-G|rystallographically equivalent 
at room temperature,ll3~E3 so that (static) charge seg- 
regation between the V atoms does not, in fact, occur. 
This result is consistent with ^^V NMR investigations 
which showed the presejrice of only one type of V atom 
at room temperature.E2l'c3 This compound is thus for- 
mally a mixed-valent dP'^ system, which has been con- 
sidered in a one-electoupbaHd picture to be a quarter- 
filled ladder compound. t3c3 We note that from mod- 
eling optical excitations in the energy range 4 meV- 
4 eV, Damascelli and coworkers initially concluded that 
the room-temperature structure of NaV2 05 is acentric;Ej 
their analysis was consistent with the V atoms on a rung 
of a ladder having oxidation states of 4.1 and 4.9, respec- 
tively. However, this group subsequently explained that 
length- and/or time-scale-of- measurement issues may be 
involved in their interpretation, such that charge dispro- 
portionation between V atoms may only occur locally 
and possibly dynamically, which could then be consistent 
with the (average long-range) cxMstal structure refine- 
ments and NMR measurements £j Theoretical suppopt 
for this scenario was provided by Nishimoto and Ohta.Ea 
Factor group analyses of the possible IR- and Raman- 
active phonon modes and comparisons with experimental 
observations at room temperature are consistent mitb-jthej 
centrosymmetric space group for the compoundH3S~E£l 
A first-principles electronic structure study based on the 
density functional method within the generalized gradi- 
ent approximation showed that the total energy of the 
centric structure is about 1.0 eV /(formula unit) lower 
than that of the acentric structure]^ consistent with the 
recent structural studies. 

One might expect that the hole-doping which occurs 
upon replacing Ca in CaV205 by Na would result in 
metallic properties for NaV205, because of the nonin- 
tegral oxidation state of the V cations and of the crystal- 
lographic equivalence of these atoms. However, NaV205 
is a semiconductor .EJ This has been explained by the for- 
mation of V-O-V molecular clusters on the rungs of the 
two-leg ladders, again resulting in a Mott-Hubbard insu- 
lator due to the on-site Coulomb repulsion ,EIrE3 where in 
this case a "site" is a V-O-V molecular cluster. Thus a 
nonintegral oxidation state and crystallographic equiva- 
lence of transition metal atoms in a compound are not 
sufficient to guarantee metallic character simply by sym- 
metry; all nearest-neighbor pairs, triplets, of transi- 
tion metal atoms must also be crystallographically equiv- 
alent, which is not the case in NaV205, since a V2 pair 
on a rung is not crystallographically equivalent to one 
on a leg in the two-leg ladders. In contrast, all V atoms 
and pairs of V atoms in mixed-valent fee LiV2 04 are 
respectively crystallographically equivalent, resulting in 
metallic character as demanded by symmetry. 

The V-O-V rung molecular clusters which are coupled 
along the ladder direction in NaV205 may be consid- 
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ered tp-|£ojjin an effective S* = 1/2 one-dimensional (ID) 
chain. Elf E3i23 Experimental support for this picture, of- 
ten quoted in the literature, is that the magnetic sus- 
ceptibility (above Tc, see below) is in agreement with 
the Bonner-Fisher prediction for the 5 p. 1/2 Heisenberg 
chain, as reported by Isobe and UedaEj Angle-resolved 
photoemission spectroscopy-^ARPES) measurements on 
NaV2 05 by Kobayashi et aZo showed that the electronic 
structure is essentially ID, despite the ostensibly 2D na- 
ture of the trellis layer, with dispersion in the oxygen and 
copper bands (below the Fermi energy) occurring only in 
the ladder direction (6-axis) . Interestingly, the dispersion 
in the lowest binding energy part of the occupied Cu Hub- 
bard band showed a lattice periodicity of 2b, which may 
reflect dynamical short-range AF and/or crystallographic 
ordering in the ladder direction. Temperature-dependent 
ARPES measurements on Nao.96V205 by the same group 
from 120 to 300 K showed evidence for the pirfidicted spin- 
charge separation in ID magnetic systems.E3 

A phase transition occurs in NaV205 at a critical tem- 
perature Tc « 33-36 K, below which the spin suscep- 
tibility —> as T ^ and a lattice distortion 
occursBj'Ea'EZi The lattice distortion results in a super- 
cell with lattice parameters 2a x 2& x 4cHj Therefore the 
transition was initially characterized as a possible spin- 
Peierls transition, which by definition is driven by mag- 
netoelastic (spin-phonon) coupling, and in which there 
is no change in the charge/spin distribution within the 
rungs/ V-O-V molecular clusters. The superstructure in 
the a and c directions, perpendicular to the V chains 
which run in the b direction, would be a result of the 
phasing of the distortions in adjacent chains/ladders. In 
this interpretation, and within the adiabatic approxima- 
tion (discussed later), one would expect that the mag- 
netic properties above Tc should be close to those of the 
5 = 1/2 Heisenberg uniform chain, and of an 5* = 1/2 
alternating-exchange Heisenberg chain below Tc. 

It has become clear, however, that the phase transi- 
tion occurring at Tc in NaV2 05 is accompanied by charge 
ordering, in contrast to a classic spin-Peierls transition. 
Therefore, the magnetoelastic coupling may only play a 
secondary role, and the spin gap may be a secondary 
order parameter. In particular, ^^V NMR experiments 
showed the presence of (inequivalent) and be- 
low— Tc, whereas only one V species was present above 
Tc.c^l This result is consistent with the solution o£-the su- 
perstructure below Tc by Liidecke and co-workersll2l using 
synchrotron x-ray diffraction. Liidecke et al. found that 
there are modulated and unmodulated chains of V atoms 
below Tc, tentatively assigned to magnetic and nonmag- 
netic chains. One interpretation of the results is that the 
^1 Y-i-4 cations segregate into alternate two-leg ladders 
which are isolated from each other within the V2O3 trellis 
layer by intervening two-leg ladders containing only non- 
magnetic V"'"^.I13 The anomalous strong increase in the 
thermal conductivity below Tc may also be due to charge 
ordering .E3 From ultrasonic measurements of shear and 
longitudinal elastic constants, Schwenk and co-workers 



have suggested that the jCtarge ordering is of the zig-zag 
type within each ladder In each of these scenarios for 
charge ordering, static charge disproportionation occurs 
such that 1/2 of the V atoms have oxidation state -1-4 
and the other half -t-5, consistent with the average for- 
mal oxidation staie of -1-4.5 in the compound. 

Koppen et aLc2l have concluded from thermal expan- 
sion measurements that the phase transition at Tc ac- 
tually consists intrinsically of two closely-spaced phase 
transitions separated by < IK, where the upper transi- 
tion is thermodynamically of second order whereas the 
lower one is first order. However, a double transition 
was not found in their specific heat measurements on the 
same crystal, which they attributed to the 50 mK tem- 
perature oscillations required by their ac measurement 
technique which were thought to broaden the two tran- 
sitions and render them indistinguishable. 

The nature of the possible charge ordering pattern has 
been studiod theoretically by several groups. Seo and 
Fukuyamaliil predicted (at T = 0) a static zig-zag chain 
of V^^ ions on each two-leg ladder, with an interpene- 
trating zig-zag chain of V"*"^ ions. They proposed that 
pairs of V"*"^ spins, one each on adjacent ladders, would 
form spin singlets, resulting in the observed spin gap. A 
similar zig-zag charge configuration_in each ladder was in- 
ferred by Mostovoy and Khomskiijlj with subsequent ex- 
perimerital support by Smirnov et a/.,E3 and by Gros and 
Valenti.O Motivated in part by the above thajmal expan- 
sion measuxement results of Koppen et aL,c3 Thalmeier 
and FuldeH proposed that the charge ordering transition 
would result in one linear chain of V"*"^ and one linear 
chain of V"^^ on each two-leg ladder, thereby then al- 
lowing a conventional spin-Peierls transition to occur at 
a slightly lower temperature, resultiiig in a double tran- 
sition as reported by Koppen et aZ.E3 A siffljlar picture 
was put forpard by Nishimoto and Ohta.E3 Thalmeier 
and Yareskoo have extensively discussed the linear-chain 
and zig-zag scenarios for charge ordering, and in addition 
have considered the alternating two-leg ladder charge 
dering pattern of the type suggested by Liidecke et al 
They point out that in both the linear and zig-zag pat- 
terns, a secondary spin-Peierls dimerization or spin ex- 
change anisotropy (in spin space) may be necessary to 
give a spin gap, whereas the two-leg ladder ordering has 
a spin gap even with no lattice distortion. Thalmeier and 
Yaresko describe the characteristic signatures of each of 
the charge-ordered models to be compared with experi- 
mental inelastic neutron scattering measurements. Riera 
and Poilblanc have discussed the influence of electron- 
phonon couplinff-pn the derived charge- and spin-order 
phase diagrams .cZl 

We have carried out x(T) measurements from 2 to 
750 K on single crystals of NaV205 along the ladder (6 
axis) direction to further characterize and clarify the na- 
ture of the magnetic interactions and ordering below and 
above Tc. We find that the magnetic properties above 
Tc are not accurately described by the 5=1/2 Heisen- 
berg uniform chain prediction with a T-independent J, 
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although a mean-field ferromagnetic interchain coupling 
can explain these data. Using our theoretical xio^^T) 
fit function for the AF alternating-exchange chain be- 
low Tc, we find that S{0) = (1 - a)/(l + a) ^ 0.034(6) 
and that the zero-temperature spin-gap of NaV2 05 is 
A(0)/fcB = 103(2) K. The 6{T) and A(r) below Tc are 
extracted. A spin pseudogap is found to occur above Tc 
with a rather large magnitude. From our specific heat 
measurements on two crystals, we find that the magnetic 
specific heat at low temperatures T < 15K is too small 
to be resolved experimentally, and that the spin entropy 
at Tc is too small to account for the entropy of the tran- 
sition. A quantitative analysis shows that at least 77% 
of the entropy change at Tc due to the transition(s) and 
associated order parameter fluctuations must arise from 
the lattice and / or charge degrees of freedom and less than 
23 % from the spin degrees of freedom. 



C. Plan of the paper 

The rest of the paper is organized as follows. Our 
notation for the Heisenberg spin Hamiltonian and for 
the reduced susceptibility, temperature and spin gap 
are given immediately in Sec. ||. Some general fea- 
tures of the high-temperature series expansion (HTSE) 
for x{T) and C(T) of S = 1/2 Heisenberg spin lattices 
and the low-temperature limits of these quantities for 
one-dimensional (ID) systems with a spin gap are then 
given. We then specialize to the S = 1/ 2 AF alternating- 



exchange Heisenberg chain in Sec. H C , where we discuss 
the HTSEs, the spin gap and the one-magnon dispersion 
relations E{A,k). In the latter subsection, we derive 
a one-parameter approximation for -©(A, k) which cor- 
rectly extrapolates to the a —^ limit and which we will 
need in order to later fit the TMRG xC^) data to extract 
A(a). We also show that the expressions for the low-T 
limits of both xC^) ^nd C{T) depend only on the spin 
gap (in addition to T). In Sec. Ill, we discuss overall 
features of the xiT) and C{T) for the uniform chain and 
then focus on the low-T behavior. The explicit forms 
of the logarithmic corrections previously found for x{T) 
are first discussed. WCpShow that a low-T expansion of 
the theory of LukyanovO gives the same first two correc- 
tions, and in addition gives the next higher-order term. 
We then compare the Bethe ansatz x{T) resultsO directly 
with the theory with no adjustable parameters or approx- 
imations. Logarithmic corrections are also found to b| 
important to accurately describe the Bethe ansatz dat. 
for C(T). We show that the lowest order correction is 
not sufficient to fit the data, and we derive the next two 
higher-order corrections by fitting the data at very low 
temperatures. 

General features of our scheme to fit numerical x(T) 
data are described in Sec. IV A, followed by a fit to the 



form chain. Due to the special requirements of, and con- 
straints on, the two-dimensional fit function necessary to 
accurately fit x(a, T) data for the alternating-exchange 
chain ove r larg e ranges of both a and T, a separate sec- 
tion. Sec. IV E| , is devoted to formulating and discussing 
this fit function. Using a fit function similar to that 
used to fit numerical x(^) data, in the next section an 
extremely accurate and precise fit is obtained over 25. 
decades in temperature to the Bethe ansatz C(T) data.El 
Our QMC and TMRG x{T) data for the alternating- 
exchange chain are presented and fitted in Sec. using 
as end-point parameters those determined for the uni- 
form chain and the dimer, respectively. The spin gap 
A{a) is extracted for 0.8 < a < 0.995 by fitting the 
TMRG_x(a, T) data at low temperatures in Sec. 0. Sec- 



tion 



VII contains a comparison of our numerical results 



with previous work. The A(a) values are compared with 
previous numerical results and with the theoretic al pre- 



diction for the asymptotic critical behavior in Sec. VII A 



Our x(T) calculations are shown in Sec. VII B to be in 
good agreement with the previous numerical results of 
Barnes and Rierafl for 0,2 < a < 0.8. The numerical 
calculations of BulaevskiiQ have been extensively used in 
the past by experimentalists to fit the xC^) of spin-Peierls 
compounds, but up to now a detailed analysis of the pre- 
dictions of this theory h as not been given. We present 
such an analysis in Sec. VII C and compare our results 
with these predictions. 

We begin the experimental part of the paper by study- 
ing the anisotropic magnetic susceptibility of a high qual- 
ity NaV205 single crystal in Sec. VIII A, where literature 
data on the anisotropy of the g factor and Van Vleck 
susceptibility are compared with our results. In the fol- 
lowing sections we illustrate the utility and application 
of many of the theoretical results derived and presented 
previously in the paper. In Sec. VIII B we present ex- 
perimental x(^) data for single c rystals of NaV205 and 
model these data in detail in Sec. VIII C using our QMC 
and TMRG x(^) data fit function for the AF alternating- 
exchange Heisenberg chain. We show that qualitatively 
and quantitatively new information about the tempera- 
ture dependences of the spin dimerization parameter and 
spin gap below Tc can be obtained from our modeling. 
This analysis also shows that spin dimerization fluctu- 
ations and a spin pseudogap are present above Tc, and 
we quantitatively determine their magnitudes. Our spe- 
cific heat measurements of NaV2 05 single crystals and 
our extensive modeling of these data are presented in 
Sec. VIII D, where we obtain quantitative limits on the 
relative contributions of the lattice, spin and charge de- 
grees of freedom to the change in the entropy due to the 
transition at Tc and to associated order parameter fluctu- 
ations. A summary and concluding discussion are given 
in Sec. fx. 



exact x(^) for antiferromagnetic Heisenberg dimer 
and two fits to the numerical xiT) dataEl for the uni- 
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II. THEORY 

In this paper we will only be concerned with the spin 
S — \/2 antiferromagnetic (AF) Heisenberg Hamiltonian 



n 



<i]> 



Jij ^ i 



(1) 



where Jij is the Heisenberg exchange interaction between 
spins Si and Sj and the sum is over unique exchange 
bonds. A > corresponds to AF coupling, whereas 
Jij < refers to ferromagnetic coupling. Note that mai 



netic nearest neighbors Sj of a given spin Si in Eq. (|1|) 
need not be crystallographic nearest neighbors. A mag- 
netic nearest neighbor of a given spin is any other spin 
with which the given spin has an exchange interaction. 

For notational convenience, we define the reduced spin 
susceptibilities x* and x* j reduced temperatures t and t 
and reduced spin gaps A* and A* as 



= ML 

Jmayi 



X 



xJ 



t = 



J 



A* 



A 



J" 



A 
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(2) 



(3) 



(4) 



where J'^'^^ and J are, respectively, the largest and aver- 
age exchange constants in the system, N is the number of 
spins, g is the spectroscopic splitting factor appropriate 
to the direction of the applied magnetic field relative to 
the crystallographic axes, and /ie is the Bohr magneton. 



A. High-temperature series expansions for the spin 
susceptibility and magnetic specific heat 

For any Heisenberg spin lattice (in any dimension) in 
which the spins are magnetically equivalent, i.e. where 
each spin has identical magnetic coordination spheres, 
the first three to four terms of the exact quantum me- 
chanical high temperature series expansion of x* (t) have 
the same form, with a particularly simple|-.forpi| if the 
series is inverted.Q For S = 1/2, one obtaingcl^E3'c3 



1 



dn 



IT. 

n=0 



(5a) 



io — 1, di — —— \ Jij, d2 — K ^ Jfi , 

4 Jmax Z-^ ■> ' g Jmax2 Z-^ n ' 

(5b) 

dz = ^ y 4 . (5c) 



Equation (^) is universal, but Eq. (^) holds only for 
spin lattices which are not geometrically frustrated for 
AF ordering and in which the magnetic and crystallo- 
graphic nearest neighbors of a given spin are the same. 
Geometrically frustrated lattices typically contain closed 
triangular exchange paths within the spin lattice struc- 
ture, such as in the 2D triangular lattice or in the 3D 
B sublattices of the fee AB2O4 oxide normal-spinel and 
A2B2O7 oxide pyrochlore structures. The uniform and 
alternating-exchange chains considered in this paper are 
not geometrically frustrated, and the magnetic and crys- 
tallographic nearest neighbors of a given spin are the 
same. It has been foundo that the terms to 0{l/t^) on 
the right-hand-side of Eq. (5a) are sufficient to quite ac- 
curately describe the susceptibilities of a variety of non- 
frustrated zero-, one-, and two-dimensional S = 1/2 AF 
Heisenberg spin lattices to surprisingly low temperatures 
t ^ 1. Higher order dn/f^ terms with n > 4 are de- 
pendent on the structure and dimensionality of the spin 
lattice. The Weiss temperature 9 in the Curie- Weiss law 
x{T) — C / {T — 9) is given by the universal expression 

9 = -diJ'^^'^/fcB. 

Because the spin susceptibility and the magnetic con- 
tribution C{T) to the specific heat can both be expressed, 
via the fluctuation-dissipation theorem and the Heisen- 
berg Hamiltonian, respectively, in terms of the spin-spin 
correlation functions, therfi_is a close relationship be- 
tween these two quantities.E3 In particular, just as there 
is a universal expression for the first three to four HTSE 
terms for x{T) of a Heisenberg spin lattice as discussed 
above, a universal expression for the first one to two 
HTSE terms for C{T) ^^^^aph a spin lattice exists and 
is given for S — 1/2 by 



NkT, 



_3_ 

32 



J2j dfj 



t^J" 



2i3jmax: 



+ 



(f) 



(6) 



The sums are over all magnetic nearest-neighbor bonds 
of any given spin Si. The first term is universal but the 
second term holds only for geometrically nonfrustrated 
spin lattices in which the crystallographic and magnetic 
nearest-neighbors of any given spin are the same. Higher 
order terms all depend on the structure and dimension- 
ality of the spin lattice. 

A common misconception is that C = if the magnetic 
susceptibility of a local-moment system obeys the Curie- 
Weiss law. This is only true classically. For Heisenberg 
spin lattices, one can easily show that the Weiss tem- 
perature 9 in the Curie- Weiss law arises from the first 
HTSE term [©(l/t)] of the magnetic nearest-neighbor 
spin-spin correlation function, which is the same, quan- 
tity that the first HTSE term of C(t) arises from.a Thus, 
e.g., for S* = 1/2 Heisenberg spin lattices at temperatures 
t ^ 1 at which the Curie- Weiss law holds, the magnetic 
specific heat is given by the universal first term of Eq. (H). 
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B. Low-temperature limit of the spin susceptibility 
and specific heat of ID systems with a spin gap 

Magnetic susceptibility. For one-dimensional (ID) 
S = 1/2 Heisenberg spin systems with a spin gap such as 
the S = 1/2 two- leg ladder (and the alternating-exchange 
chain), Troyer, Tsunetsugu, and WiirtzEil derived a gen- 
eral expression for x* (t) which approximately takes into 
account kinematic magnon interactions, given by 



^^^'^~t 1 + Mt) 



Jo 



d{ka) 



(7a) 



(7b) 



where Sk = E{k)/J"^^^, E{k) is the nondegenerate one- 
magnon (triplet) dispersion relation (the Zeeman degen- 
eracy is already accounted for) and a is the (average) 
distance between spins. This expression is exact in both 
the low- and high-temperature limits. For the isolated 
dimer, for which = A* = 1, Eq. ( [7a| ) is exact at all 
temperatures. Inserting z{t) = e~^/* for the dimer into 
Eq. ( [7a| ) yields the correct result 



t l + 3c-i/t ' 



(dimer) 



(8a) 



(8b) 



The x*(^) in Eq. ( pa]) for the antiferromagnetic Heisen- 
berg dimer is plotted in Fig. Q the fit shown in the figure 



will be presented and discussed later in Sec. |IVB . 

At low temperatures t <^ A* and t ^ one-magnon 
bandwidth/ J™^'^, and for a dispersion relation with a 
parabolic dependence on wave vector k near the band 
minimum 



m 

jmax 



'(kaf 



(9) 



one can replace Sk in Eq. (7b) by the approximation (|^) 
and replace the upper limit of the integral in Eq. (7b) by 
oo, yielding z{t) — e~^*/*\/t/(2\/7rc*). Substituting this 
result into Eq. ( [7a| ) gives the low-t liniitEII 



A = 



A 
In 



^=2 



(10a) 



(10b) 



This result is correct for any ID = 1/2 Heisenberg spin 
system with a spin gap and with a nondegenerate (ex- 
cluding Zeeman degeneracy) lowest-lying excited triplet 
magnon band which is parabolic at the band minimum. 
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FIG. 1. (a) Magnetic susceptibility x (°) versus tempera- 
ture T for the spin 5=1/2 Heisenberg dimer with antifer- 
romagnetic exchange constant J. The fit from Sec. [V B is 



shown by the solid curve, (b) Semilog plot of the fit deviation 
vs T. The lines connecting the points in (b) are guides to the 
eye. 

On the other hand, the low-temperature limit of x*(0 
for the isolated dimer in Eq. ( pq ) is of the same form 
as Eq. (10a), but with 7 = 1. Thus, for ID systems 



consisting essentially of dimers which are weakly coupled 
to each other, a crossover from 7=1 to 7 = 1/2 is 
expected with decreasing t. 

The parameters A and 7 can be determined if very 
accurate x*(^) an d A * data are available. Taking the 
logarithm of Eq. (|lOa|) yields the low-i prediction 



ln[x*(t)] + ^ = lnA-7lnt , 



(11a) 



so plotting the left-hand-side vs lat allows these two 
parameters to be determined. Alternatively, assuming 
7 = 1/2 , one can obtain estimates of A and A* using 
Eq. (10a), according to 



Hx*Vt) = 



A* 



and / or 



dHx*Vi) 
d{l/t) 



= A*. 



(lib) 



(He) 
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Specific heat. The \ow-t limit of the magnetic coatri- 
bution C(T) to the specific heat for the same modelEJ is 
calculated to be 



C{t 



NL 



0) _ 3 / A* 
I ~ 2 

1 - 



1/2 
t 



3/2 
t 



-A'/t 



(12) 



Note that, in addition to the ratio t/A* ~ fceT/A 
of the thermal energy to the spin gap, the magnitude 
of X* in Eqs. ( [lO| ) is determined by the actual value of 
the curvature c* at the triplet one-magnon band mini- 
mum, whereas the magnitude of C in Eq. ( [T^ ) depends 
only on the ratio of c* to A*. These formulas have been 
applied in the literature to model experimental data for 
alternating-exchange chain and two-leg-spin ladder com- 
pounds. However, with one exceptionE3 to our knowl- 
edge, these modeling studies have not recognized that the 
prefactor parameter and the spin gap are not indepen- 
dently adjustable parameters. For a given spin lattice, 
they are in fact uniquely related to each other. Their re- 
lationship for the 5* = 1/2 two- leg Heisenberg ladder was 
studied in Ref. 52. For the alternating-exchange chain, 
we estimate the relationship between c* and A* below in 



Sec. [IC3 



C. Alternating-exchange chain 

The S — 1/2 AF alternating-exchange Heisenberg 
c^^^n Hamiltonian is written in three equivalent ways 

H = JlS2i-l ■ S2i + ^2'5'2i • '5'2i-|-l (13a) 
i 

= '^^JiS2i-i ■ S2i + aJiS2i ■ S2i+i (13b) 

i 

= ^ J(l + 5)S2^-l ■ S2; + J(l - 5)S2^ ' S2^+l , (13c) 

i 

where 



Ji ^ J(l + (5) = 



2J 
1 -h a 



1-5 



a = — = 



J 



J I 1 + 5' 
'II - I- -^1 ~ -^2 _ 1 - g 

J ^ 2 J ~ l + a 
J1 + J2 l + a 

7i — Jl — 7^ — I 



(14a) 
(14b) 
(14c) 
(14d) 



with AF couplings Ji > J2 > 0, < (a, (5) < 1. 
The uniform undimerized chain corresponds to a 
1, 5 = 0, J i = J2 = J. The form of the Hamiltonian 
in Eq. (13c) is most appropriate for chains showing a 



second-order dimerization transition at Tc with decreas- 
ing T. If the exchange modulation (5 <C 1 (a ~ 1), the 
(average) J below Tc can be identified with the exchange 
coupling in the high-T undimerized state. 

In spin-Peierls systems, the spin-phonon interaction 
causes a lattice dimerization to occur below the spin- 
Peierls transition temperature, together with a spin-gap 
due to the formation of spin singlets on the dimers. 
The Hamiltonian can be mapped onto the spin Hamil- 
tonian (13) (with renormalized exchange constants) only 
in the adiabatic parameter regime, in which the relevant 
phonon energy is much smaller than J. If such a map- 
ping cannot be made, dynamical phonon effects (quan- 
tum mechanical fluctuations) become important and the 
x(r) can be significanjtiy .-eLifferent from that predicted 
from Hamiltonian (^3|).eJ"E3 This issue will be discussed 
furt her wh en modeling the x{T) data for NaV205 in 



Sec. VIHB 



1. High-temperature series expansions 



Magnetic Susceptibility. For the alternating-exchange 
chain, according to our definition one has J™^'^ = J^. 



Then using the definition for a in Eq. (14b), the d„ HTSE 
coefficients in Eqs. ( ^b| ) and ( |5^ ) become 



do = 1, di ^ — - — , d2 
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,^3 = 4^. (15) 



One can change variables from a and Ji in x*(a, t) to 5 
and J in x*i5,t) using Eqs. (nj) which give 



t = 



1 + 5 ' 



(16a) 



(16b) 



We write the resulting HTSE for the inverse of x* {5, t) 
as 



1 



4x* t 



n=0 



dn 



(17a) 



where we find 



-J 1 ^ 1-3 1 + *^^ -7 l + 3<5^ 

do = 1, di = -, da = — - — , da = — — — . (17b) 

An important feature of this HTSE of x*{^:'t) is that it 
is an even (analytic) function of 5 for any finite temper- 
ature. This constraint must be true in general and not 
just for the terms listedjj because x*((5, t) cannot depend 
on t he si gn of i5 = (Ji — J2)/(2J): the Hamiltonian in 
Eq. (13c) is invariant upon such a sign change. Physi- 
cally, a negative 5 would simply correspond to relabeling 
all Si — > Si+i, which cannot change the physical prop- 
erties. We will use this constraint that x*((5, t) must be 
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an even function of S to help formulate our fitting func- 
tion (after a change back in variables) for our QMC and 
TMRG x*{'^it) calculations for the alternating-exchange 
chain. This constraint is important because it allows a 
fit function for x*(a,i) to be formulated which is accu- 
rate for a 1 ((5 ^ 1), a parameter regime relevant 
to compounds exhibiting second-order spin-dimerization 
transitions with decreasing temperature. 

Magnetic specific heat. Using J"^'^^ = and a = 
J2/J1, the general HTSE expression in Eq. (||) yields 
the two lowest-order HTSE terms for the magnetic spe- 
cific heat C{T) of the S* = 1/2 AF alternating-exchange 
Heisenberg chain as 



Dispersion Relations for Alternating Ciiain 



C[t) 

NkB 



3 rl 
32 .~ 



a 



2t^ 



O 



(?) 



(18) 



2. Spin , 



The spin gap A* (a) of the alternating-exchange chain 
was determined to high (< 1%) accuracy for < a < 0.9, 
in a increments of 0.1, using multipprision methods by 
Barnes, Riera, and Tennant (BRT).e3 They found that 
their calculations could be parametrized well by 



A*(a) 



A(«) 
Ji 



A*{S) = 



(l-a)3/4(l + a) 



1/4 



253/4 



(19a) 



(19b) 



The same A*G£) was found in numerical calculations 
by Ladavac et al!t3 for 0.01 < S < l,jchereas calculations 
for 0.03 < S < 0.06 by Augier et alsB yielded som ewhat 
smaller values of A* than predicted by Eq. ( 19b ). The 
asymptotic critical behavior of A* as the unifpxipJimit 
is approached (a — *■ 1, 5^0) has been givenOoiEl as 



A*{S) cx 



ln(5|i/2 



(20) 



thus the parametrization in Eq. (19b) evidently indicates 
that the fitted data do not reside within the asymp- 
totic critixial regime. Alternatively, Barnes, Riera, and 
TennantH suggested that Eq. (20) may not be the cor- 
rect form for the asymptotic critical behavior. On the 
other hand, Uhrig et al£3 fitted their T = density 
matrix renormalization group (DMRG) calculations of 
A*(5) for 0.004 < 6 < 0.1 to a power-law behavior with- 
out the log correction and obtained A* « 1.57 (5"'^^. We 
will further discuss the above spin gap calculation results 
later in Sec. VII A after deriving our own A* (a) values 
from our TMRG x*{a,t) data in Sec. ^. 



E(a,k)/J = [A*(a)2 + f{a) sin2{ka)]^'2 




0.4 0.6 

ka/7t 

FIG. 2. Dispersion relations E{a,k) for the S = 1/2 an- 
tiferromagnetic alternating-exchange Heisenberg chain. The 
dashed curves for alternation parameters a = 0, 0.2, 0.4, 0.6 
and 0.8 are the dimer series expansion results of Barnes, Ri- 
era and Tennant (Ref. ^3|) , the solid curves for these a values 
are from our expression given in the figure and in Eqs. (p3[), 
and the solid curve for a = 1 is the known exact result for 
the uniform chain, which by construction is the same for the 
expression in the figure at this a value. 



3. One-magnon dispersion relations 

Barnes, Riera, and Tennant have computed the dimer 
series expansion of the dispersion relation e{a, k) = 
E{a,k)/Ji for the one-magnon {S = 1) energy E{a,k) 
vs wave vector k along the chain for thcJowest-lying one- 
magnon band, up to fifth order in q;,o which we write 
as 



e(a. A;) = ^^a„(a) cos(2nfco) , 

n=0 



(21) 



where a is the (average) spin-spin distance, which is 1/2 
the lattice repeat distance along the chain in the dimer- 
ized state. Plots of e(a, k) for a — 0, 0.2, 0.4, 0.6, and 
0.8 up to fifth order in a, as given in Fig. 4 of Ref. |5^, 
are shown as the dashed curves in Fig. ^j. The curves 
are symmetric about ka = 7r/2, so the same spin gap 
A* (a) = X]^o '^nict) occurs at fca = and tt. This fifth- 
order approximation yields A* (a) values for a < 0.9 in 
rather close agreement with BRTs' results discussed in 
the previous section. For a dimer series expansion we ex- 
pect the average energy of the one-magnon band states 
to be nearly independent of a, i.e.. 



1 



e{a, ka) d{ka) 



1 



(22) 



Indeed, upon inserting BRTs' fifth order expansion coef- 
ficients into Eq. ( |2l| ) and the result into Eq. ( p2|) , we find 
that this sum rule is satisfied to within 1% for < a < 1. 

Also shown as a solid curve in Fig. ^ is the exact resiilt 
e{k) = (7r/2)| sin(fca)| for the uniform chain (a = 
This e{k) has a cusp with infinite curvature (at ka = Q 
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FIG. 3. The function /(A*) in the expression in the figure 
and in Eqs. ( p^ ) for the one-magnon dispersion relation of the 
S = 1/2 antiferromagnetic alternating-exchange Heisenberg 
chain, where A* = A/Ji is the spin gap. The dependence 
/(a) is also plotted for the assumed form of A*(q) shown. 

and tt) which cannot be accurately approximated by a 
Fourier series with a small number of terms. This sin- 
gular behavior is evidently closely related to the critical 
behavior of A* (a 1) discussed above. In order to later 
model our TMRG x*(a,i) data close to, but not in, the 
low-t limit, we will need an expression for £(a, k) which 
is correct in the limit a — *■ 1 and which also reproduces 
reasonably well the e(a, k) of BRT. We found that the 
simple one-parameter (A* ) form suggested earlier hy one 
of us in the context of the S = \/2 two- leg laddero 



e{A*,k) = [A*^ + /2(A*) sin2(fca)]i/2 



(23a) 



is satisfactory in these regards for the AF alternating- 
exchange chain over the entire range < A* < 1. The 



function /(A*) is determined here by the sum rule (22) 
which yields the condition 



E 



/2(A*) 



2A* 



(23b) 



where E(a;) is the co mple te elliptic integral of the sec- 
ond kind. From Eq. ( |23b| ), / varies nonlinearly with A* 
from /(A* = 0) = it/2 to /(A* = 1) = 0, as shown in 
Fig. 1^. From an ind epen dently determined dependence of 
A* on a as in Eq. (|l9a[) , one can then determine /(a) as 
also shown in Fig. ^ Using the fifth-order A* (a) values 
of BRT in Fig. the resulting dispersion relations (|2^ ) 
for a = 0, 0.2, 0.4, 0.6, and 0.8 were calculated and are 
shown as the solid curves in Fig. H, where they are seen to 
be in close agreement with the respective dashed curves 
of BRT. An important difference for large a, however, is 
that the £(A*,fc) in Eqs. ( p3| ) properly reduces by con- 
struction to the exact e(a, fc) for a 1, whereas the one 
in Eq. ( pi]) with a finite number of terms does not. 

Close to the one-magnon band minimum, the square 
root and the sine function in the dispersion relation in 
Eq. (23a) can be expanded, yielding 



e(A*,fc-> 0) w A* 



1 /'(A* 
2 



A* 



-{ka) 



(24) 



A comparison of Eqs. (2^) and 
rameter c* in the formu 



) shows that the pa- 
as for x*(i ^ 0) [Eqs. (0)] and 
C(t ^ 0) [Eq. (p^)] is a unique function of A* which in 
our approximation is given by 

1 /^(A* 

2 A* 

with /2(A*) given by Eq. (^3b| ). Thus both x*(t 0) 
and C{t 0) for the alternating-exchange chain only 
depend on the single parameter A* (in addition to t). 
Explicitly, we obtain 

1 /A*\i/2 



^(A*) 



(25) 



(t) 



-A*/t 



(26) 



/27r/(A*) 

As might have been anticipated, the only thermodynamic 
variable is the ratio t/A* = fceT/A of the thermal en- 
ergy to the spin gap. The numerical prefactor depends 
explicitly (only) on the reduced spin gap A* = A/Ji. 
Similarly, the magnetic specific heat is obtained as 

C{t^O) 3 A* /A*n3/2 



2^ /(A* 



t 



3/ t 
iVA^ 



(27) 



where again the same characteristics are present as just 
discussed for x* (t) ■ The variations of the prefactors with 
A* for x*{t) ^^^d C{t)ca.n both be ascertained from the 
plot of /(A*) in Fig. |3[ In particular, when a ^ 1 (5 ^ 
1) for which A* ^ 1, / is nearly a constant. For our 
and our readers' convenience when modeling materials 
showing small spin gaps, we have fitted our numerical 
/(A*) calculations for < A* < 0.4 by a third-order 
polynomial to within 2 parts in 10*, given by 

/(A*) = I - 0.034289A* - 1.18953A*^ + 0.40030A*^ 

(28) 

By a chaiige in variables to (J, 5) and using the A* (a) 
in Eq. (^9|)£j we obtain the following forms which are 
more useful for modeling materials with small spin gaps, 
especially those showing second-order spin dimerization 

transitions with decreasing T: 

1 /A*\i/2 
r{t - 0) = ' f ^) e-^-/* . (29a) 



(l + <5)x/27r/(A*) V t 



C{t ^ 0) 
Nkn 



A* /A*\3/2 

(l + 5)\/2^7W)' 



1 J_ -(—Y 

A* 4 V A^/ 



(29b) 



/(A*) = I - 0.033933A* - 1.19607A*^ -I- 0.92430A*^. 

(29c) 

Note that in these formulas A* /t = A*/t = A/{kjiT). 
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III. THEORY: S 



1/2 UNIFORM HEISENBERG 
CHAIN 



A. Magnetic spin susceptibility 

The uniform S ~ 1/2 chain is one Umit of the 
alternating-exchange chain with Jij = J, a = 1, (5 = 0, 
and with no spin gap [the x*{t 0) and C{t ^ 0)/t are 
finite]. The spin susceptibihty was calculated accurately 
by Eggert, Affleck and Takahashi in 1994,0 and recently 
refined by KliimperB as shown in Fig. ^(a) where only the 
calculations up to t = 2 arc shown. An expanded plot of 
the data for t < 0.02, including the exact value I/tt^ at 
i = 0, is shown in Fig. ^(b), along with a fit (Fit 2) to 



rV^ . The 



the data to be derived and discussed in Sec. 
most recent calculations of Ref. ^ have an absolute ac- 
curacy estimated to be w 1 x 10^^ and show a broad 
maximum at a temperature T"^^^ with a value By 
fitting data points near the maximum by up to 6th order 
polynomials, we determined these numerical values to be 
given by 



j,max ^ 0.6 408 510(4) J/ZCB 



X 



^maxj 
WJ4 



max^^max 



0.146 926 279(1) , 



0.0 941 579(1) 



(30a) 



(30b) 



(30c) 



These values are consistent within the errors with those 
found by Eggert et al.p but are much more accurate. For 
one mo le of spins, setting N = Nj^ (Avogadro's number) 
in Eq. (^0^) yields 



X 



maxymax ^ 0.0 353 229(3); 



cm^ K 
mol 



(31) 



Note that the product ;^max2.max (|30cD and (|3|) 

is independent of J, and hence is a good initial test of 
whether the 5" = 1/2 AF uniform Heisenberg chain model 
might be applicable to a particular compound. 



1. High-temperature series expansions 
The coefficients c„ of the HTSE for x*(t), 

OC 

n=0 

are given up to 0{l/t'^) 



1 1 

Co = 1 , Cl = -- , C2 = , C3 = — , C4 



5 

384 



7 



C5 = 



1280 



C6 = 



133 
30720 



1 



4032 



(32b) 
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FIG. 4. (a) Magnetic susceptibility x versus temperature 
T for the spin S = 1/2 nearest-neighbor antiferromagnetic 
Heisenberg chain (Ref. ^) (•). (b) Expanded plot of the data 
in (a) for < t < 0.02 , together with a fit (Fit 2, solid curve) 
obtained in Sec. IVC to the data of Kliimper and Johnston 
(Ref The fit is not shown in (a) because on the scale of 
this figure the fit is indistinguishable from the data. 

Inverting the series, we obtain the corresponding rf„ co- 
efficients in Eq. (|5^) as 



do = 1 , rfi = i , d2 
29 



4' ^^-12' ^^-128' 



3840 



317 
92160 



d7 



11 



71680 



(33) 



The dn coefficients with n = 0, 1, 2, and 3 arc of course 
in agreement with Eq. (|l^) for a = 1. 



2. Logarithmic corrections at low temperatures 

At low temperatures a simple expansion of thermo- 
dynamic properties in for instance the variable t is not 
possible. Such a nonanalyticity in t can be viewed as 
due to the strong correlation of the quasiparticles, i.e., 
the elementary excitations of the system are not strictly 
free; they show rather nontrivial scattering processes. 
Spinous with low energies ei and £2 have a scattering 
phase 0(ei, £2)^ 0o + const/| log(eie2)|. From this prop- 
erty it is cleaiEa that an expansion in the single variable 
t is not possible, but has to be supplemented by a term 
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l/log(t). Although being very intuitive, this physical 
picture on the basis of scattering processes of spinous 
has not played any important role in the investigation of 
logarithmic corrections until recently.d 

Logarithmic dependencies of physical quantities have 
been known for the spin- 1/2 Hcisenberg chain for a rather 
long time. Usually, a quantum chain with critical cou- 
plings leads to critical correlations only in the thermo- 
dynamic limit 1/i = and at T = 0, where L is the 
length of the chain. If one of these conditions is not met 
the physical properties receive nonanalytic contributions 
in terms of 1/L or T. From the renormalization group 
point of view the existence of logarithmic corrections is 
reflected by the perturbation of the (critical) fixed point 
Hamiltonian by some marginal operator. Such operators 
usually exist only for isotropic systems. 

The investigation of the size dependence of energy lev- 
els of critical quantum chains was started more than 
a decade ago. For the isotropic spin-1/2 Heisenberg 
chain, expansions in 1/L and additional logarithmic cor- 
rections (l/LlogijCtc.) were found in lattice approaches 
(Bethe ansatzt3~l£a) as well as in field theory [RG study 
of the Wess-Zumino-Witten (WZW) model with topolog- 
ical term k = 1 (Refs. |l|,|9|,^]. 

Many of these earlier results are still relevant for the 
issues discussed in this section due to an equivalence of 
many-particle systems at T = 0, 1/L > (groundstate 
properties of finite chains) and those atT>0, 1/ L = {) 
(thermodynamics of the bulk). This leads to asymp- 
totic series where T and 1/L play very similar roles. To 
our knowledge the first explicit report on log(T) correc- 
tions in the magnetic susceptibility resulting in an infinite 
slope at T = was given in Ref. ^. Including highcr^ji^CT 
terms, the asymptotic expansion Xiti^) for X*{t) is™ 



1 

1 



In (C 



^+2Z- 



{2CY 
ln£ 



{2Cf {2Cf 



C = \n{to/t) 



(34a) 
(34b) 

(34c) 



where to is a nonuniversal (undetermined) parameter. In 
Ref. ^ the field theoretical prediction on the basis of the 
WZW model was compared with the results of thermody- 
namic Bethe ansatz calculations and showed convincing 
agreement in an intermediate temperature regime. Using 
up to the first logarithmic correction term in Eq.J34a), 
Eggert, Affleck, and Takahashi estimated « 7.7jj so at 
low temperatures t < 0.01 the parameter £ ^ 1. 

A general feature of field theoretical and lattice ap- 
proaches is their restriction to "low" and "high" temper- 
atures, respectively. Field theoretical studies suffer at 
high temperatures from the neglect of (infinitely many) 
irrelevant operators. Lattice studies show convergence 



problems at low temperatures as increasingly larger sys- 
tems have to be studied in order to avoid finite-size ef- 
fects. In addition, the comparison of field theory and 
lattice results can only verify or falsify the universal as- 
pects of an asymptotic expansion. Non-universal quanti- 
ties like to which derive from some coupling constant of 
a marginal or irrelevant operator (undetermined within 
the field theory) can at best be fitted as done in Ref. ||. 

The latter problem of determining the coupling 
constants pin an effective field theory was solved by 
Lukyanovfl who used a bosonic representation of the 
Heisenberg chain. The coupling constants were fixed by 
a comparison of the susceptibility data x{T = 0, h) ob- 
tained by him with Bethe ansatz calculations for mag- 
netic field /i at T = 0. Eventually, the x(T > 0, /i = 0) 
data could be calculated without any need of a fit pa- 
rameter. „ 

LukyanovD obtained the following analytical low- 
temperature expansion of X*(t)j 



IT 



0{9)] 



where g obeys the transcendental equation 



\\ng +- = C 

2 g 



or equivalently 



V~9 



.1/9 



to 

t 



with a unique value of to given by 



to 



e7+(i/4) ^ 2.866 257 058 



(35a) 



(35b) 



(35c) 



(35d) 



where 7 « 0.577 215 665 is Euler's constant. Lukyanov 
showed that his Xj*^. g(t) is in agreement with the numer- 
ical data of Eggert, Affleck and Takahashii at low tem- 
peratures t > 0.003. 

In the following, we will compare high-accuracy nu- 
merical Bethe ansatz calculations carried out to much 
lower temperatures by Kliimper and JohnstorO with this 
theoryQ in some detail because this theory is exact at 
low temperatures with no adjustable parameters. The 
calculations of Ref. ^ are based on lattice studies, how- 
ever without suffering from the usual shortcomings. By 
means of a lattice path integral representation of the fi- 
nite temperature Heisenberg chain and the formulation 
of a suitable quantum transfer matrix (both quite analo- 
gous to the numerical TMRG calculations presented later 
in this paper) a set of numerically well-posed expressions 
for the free energy was derived. In more physical terms 
the method can be understood as an application of the 
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familiar though often rather vague concept of quasipar- 
ticles to a quantitative description of the many particle 
system valid for all temperatures T and magnetic field 
values hu The work can be understood as an evaluatioa 
of the full scattering theory of spino us and antispinons.u 
Our iterative solution of Eq. ( p5b| ) yields the expansion 
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1 



ln£ 
'2C 



(In/:) 



ln£ 



{2Cf 



O 



1 



{2Cf 



(36) 



A lo g-log plot of 5 vs t obtained by numerically solving 
Eq. (35c) is shown in Fig. || (solid curve), along with 
its lowest-order approximation g{t) « 1/C = l/\n{to/t) 
(dashed curve). This approximation is 1.1% larger than 
the exact result at t = 10^^°, with the discrepancy in- 
creasing steadily to 5.8% at t ~ 10^^^ and 8.5% at 
Substituting Eq. 



t = 10- 



xrt,iog(^) - 



( p6D into (35a) gives 



1 

(ln/:)2 -ln/:+ (3/4) 



{2Cf 



o 



1 



V3 



1 



{2CY 
(37) 



The first three terms are identical with those in Eq. (34b), 
but the constant term in the numerator of the fourth term 



is not the same as in Eq. (34b), indicating that Eq. (34t) 
is not accurate to order higher than 0[l/(2£)^]. 

An important issue is the accuracy to which the log 
expansion approximation X\t log (^) ^1- ( |37|) reproduces 
the X\t g(*) prediction of the original Eqs. (|35|). We have 
calculated both quantities to high accuracy and plot the 
difference vs t, for the range lO"'^" < t < 0.5, in Fig. |^. 
The X\t iog(0 is seen to increasingly diverge from Xit g(t) 
with increasing t. 

When comparing the predictions of Lukyanov's theory 
with numerical results such as obtained from the Bethe 
ansatz, it is important to know at what temperature the 
low temperature expansion in Eqs. (|3^) ceases to be ac- 
curate ( "accurate" must be defined) with increasing tem- 
perature. There are three aspects of this issue that need 
to be addressed. The first and second aspects concern 
the temperatu res a t which the unknown 0{g'^) and 0{g) 
terms in Eq. (35a) become significant, respectively; we 
will return to these two issues shortly. The third aspect 
is whether the log expansion approximation Xit log(^) 
Eq. (^ can be used in this comparison. The absolute 
accuracy of the most recent Bethe ansatz calculations^ 
is estimated to be ~ 1 x 10^". From Fig. |[ we see that 
Xit iog(0 approximates x*t g(0 to this degree only for tem- 
peratures t < lO"'^" [we infer that the previous Eqs. (^) 
only apply to this accuracy at similarly very low temper- 
atures]. Therefore, to avoid this unnecessary approxima- 
tion as a source of error at higher temperatures, we will 
henceforth compare the numerical Bethe ansatz calcu- 
lations with Xit g(^) calculated from Lukyanov's original 
Eqs. (H). 
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FIG. 5. Log-log plot vs temperature T of the func- 
tion g{T) (solid curve) and its lowest-order approximation 
g{T) = l/ln(To/r) (dashed curve) in Lukyanov's theory 
(Ref. ^ for the S = 1/2 antiferromagnetic uniform Heisenberg 
chain over the temperature range 10"^" < t < 0.5. 
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FIG. 6. Log-log plot vs temperature T of the difference 
between our approximate logarithmic expansion Xi* log of 
Lukyanov's theory (Ref. ^ and his exact prediction Xit.g for 
the low-temperature limit of the magnetic susceptibility of the 
spin 5 = 1/2 antiferromagnetic uniform Heisenberg chain, for 
the temperature range 10""^" < t < 0.5. 



A comparison jaf the low-temperature Bethe ansatz 
X*it) calculationsB and Lukyanov's theoretical Xit g P^'^" 
diction is shown in Fig. 0(a). On the scale of this figure, 
the two results are identical. The (small) quantitative 
differences between them are shown as the filled circles 
in Fig. 0(b). The lower error bar on each data point in 
Fig. M(b) is 1 X 10"^ to indicate the scale. The upper error 
bar is the estimated uncertainty in Xit g arising from the 
presenc e of the unknown 0{g^) and higher-order terms 
in Eq. (35a), which was set to g^{t)/Tr'^; the uncertainty 
in the contribution, ^ \/3t^g{t)/TT^, is negligible at low 
t compared to this. At the lower temperatures, the data 
agree extremely well with the prediction of Lukyanov's 
theory. At the highest temperatures t > 10"'^, higher or- 
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FIG. 7. Semilog plots vs temperature T at low T of (a) 
numerical Bethe ansatz magnetic susceptibility (x*) data for 
the S — 1/2 uniform Heisenberg chain (Ref. |5|) (•) and the 
prediction Xit,g (solid curve) of Lukyanov's theory (Ref. 
and (b) the difference between these two results (•). In (b), 
the upper error bar is the estimated uncertainty in Xi't.g (see 
text). 

der terms also become important, as inferred from our 
empirical fits (Fits 1 and 2) below to the numerical data. 

Irrespective of the uncertainties in the theoretical pre- 
diction at high temperatures just discussed, we can safely 
conclude cirectly from Fig. ^(b) that the Bethe ansatz 
X*{t) data J are in agreement with the exact theory of 
Lukyanovfl to within an absolute accuracy of 1 x 10~^ 
(relative accuracy ~ 10 ppm) over a temperature range 
spanning 18 orders of magnitude from t — 5 x 10~^^ to 
t = 5 X 10^^. The agreement is much better than this at 
the lower temperatures. 



B. Magnetic specific heat 

The magnetic specific heat C of the S = 1/2 AF uni- 
form Heisenberg chain was recently calculated to high ac- 
curacy by Kliimper and Johnstcm over the temperature 
range 5 x 10~^^ < k-aT/J < 5B The accuracy is esti- 
mated to be 3 X 10-'^°C{t). The results for T < 2J/fcB 
are shown in Fig. |^(a). The initial T dependence is ap- 
proximately (see below) linear, and is given exactly in 
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FIG. 8. (a) Specific heat C vs temperature T (•) for the 
S — 1/2 antiferromagnetic uniform Heisenberg chain (Ref. |^. 
(b) Specific heat coefficient C/T vs T from the data in (a). 
The area under the curve in (b) from T = to T = 5J/fcB is 
99.4% of In 2. 



the t — Q limit by 



C{t 0) 

iVfcR 



2 

3^ 



(38) 



The data show a maximum with a value C"^^^ at a tem- 
perature Tq^^. By fitting 3-7 data points in the vicinity 
of the maximum by up to 6th order polynomials, these 
values were found to be 



J 

^max 



0.48 028 487(1) , 
0.3 497121235(2) 



(39) 



The electronic specific heat coefficient C{T)/T is plot- 
ted vs temperature in Fig. ^(b). As expected from 
Eq. (^), the data approach the value (2/3)iVfc|/J for 
t — > 0. The initial deviation from this constant value 
is positive and approximately (see below) quadratic in 
t. The data exhibit a smooth maximum with a value 
{C/T)'^^^ at a temperature T^^^, values which we deter- 
mined by fitting polynomials to the data in the vicinity 
of the peak to be 
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FIG. 9. Entropy S vs temperature T for the 5" = 1/2 
antiferromagnetic uniform Heisenberg chain, obtained from 
the data in Fig. ^(b). The entropy is normalized by 
S{T = oo) = NkB In 2. 



HBJ^ C/T 
J 

{C/T)""'''J 



0.30 716 996(2) 



= 0.8 973 651576(5) . 



(40) 



The magnetic entropy S(T) is determined by integrat- 
ing the C{T)/T data in Fig. |(b) vs T and the result, 
normalized hy S{T ^ oo) — NkB\n2, is plotted vs T 
in Fig. ^. This figure allows one to estimate the max- 
imum magnetic entropy that can be associated with a 
dimerization transition or any other magnetic transition 
involving S = 1/2 Heisenberg chains which are weakly 
coupled to each other [assuming that the (average) J 
does not change at the transition]. For example, for 
NaV205 where ksTc/J « 0.057, one can estimate from 
Fig. ^ that the magnetic entropy at Tc cannot exceed 
0.056i?ln2 = 0.32J/molK, where R is the molar gas 
constant. The reason this value is the upper limit is that 
magnetic critical fluctuations will increase the specific 
heat, and hence the entropy, above Tc and thus reduce 
it at (and below) Tc, by conservation of magnetic en- 
tropy, compared to the values for the isolated chain at the 
same reduced temperatures. Similarly, the C(T) data in 
Fig. |^(a) allow one to estimate the minimum lattice spe- 
ciflc heat contribution C'''*(T) above Tc if the C^^\T) has 
not been determined previously from experiments and/or 
theory directly. 

At low temperatures, the electronic specific heat co- 
efficient C(T)/T becomes independent of temperature 
(apart from logarithmic corrections, see below), as does 
the spin susceptibility x*(^)j just as in a metal (Fermi liq- 
uid). Therefore it is of interest to compute a normalized 
ratio of these two quantities. For a metal, the relevant 
quantity is the Wilson-Sommerfeld ratio, which for 
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FIG. 10. The Wilson-Sommerfeld ratio Rw = 47r^x*(*)i/ 
[3C(f)] between the spin susceptibility and electronic specific 
heat coefficient for the 5* = 1/2 AF Heisenberg chain vs re- 
duced temperature t = ksT/J. In the Wilson ratio, we have 
set kB = 1. 

5=1/2 quasiparticles reads, in the notation of this paper 
and with set to 1, 



Rw{t) 



4^'X*(0^ 
3C{t) 



(41) 



In a degenerate free electron gas, R-^ — 1 and is indepen- 
dent of t. For exchange-enhanced metals 1 < i?w ^ 10, 
for S" = 1/2 Kondo impurities in a metal the Wilson 
ratio associated with the impuritiespis i?w — 2, and for 
many heavy fermion metals Rw ~ 2Jl3 Plotted in Fig. |l^ 
is Rw{t) for the S' = 1/2 AF Heisenberg chain, where 
C[t)/t and x*[t) were given above in Figs, ^(b) and ||, 
respectively. For i ^ 0, the Wilson ratio for the = 1/2 
Heisenberg chain is exactly 2. With increasing Rw 
seen to be nearly independent of t to within ±10% up 
to t w 0.4, but the influence of the logarithmic correc- 
tions to both x(T) and C(T) are quantitatively impor- 
tant. Although the logarithmic corrections for x{T) and 
C(T) oppose each other in their ratio in Rw{t), the log- 
arithmic corrections for x(^) win out, giving a net ~ 10% 
increase in Rwit) with increasing t at low t. At higher 
t, the system crosses over to the expected local moment 
Heisenberg behavior where Rw oc t^. Thus as far as 
the thermodynamics is concerned, the uniform Heisen- 
berg chain behaves at low temperatures as expected for 
a Fermi liquid, apart from the influence of the logarith- 
mic corrections. This quasi-Fermi liquid behavior arises 
because the elementary excitations at low temperatures 
are S = 1/2 spinous which are fermions with a Fermi 
surface (i.e., Fermi points in one dimension). Since the 
spinous carry no charge, the chain is an insulator. The 
deviation of the Wilson ratio from unity and the loga- 
rithmic corrections are due to spinon interactions. 
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1. High-temperature series expansions 

The HTSE for the-specific heat of a spin S AF uniform 
Heisenberg chain isE3 



Logarithmic Specific Heat Corrections 



C{T) _ x^r 



^ Cn[X) 
/ ^ fn 



a; = 5(5+1), t 



J ' 



(42a) 



(42b) 



C4 



5040 



ci = i , C2 = ^(3- Sx- 3x2) , 



C3 = ^(3-16a;-4a;2) 



(192 - 14322; + 11232;^ + SOOx^ + leOa;"*) , 



C5 



21600 



(414 - 3768x + 6635x2 _^ 2624^^ + 480x'*) . 

(42c) 

Speciahzing Eqs. (|4|) to 5 = 1/2 (x = 3/4) then gives 
C(r) 3 



NkB 16*2 



1+y ^ 



C2 



^^ = -16 



C4 



7 

256 



C5 



(43a) 



917 
7680 ' 
(43b) 



The two C(T) HTSE terms of order 1/^2 and in 
Eqs. ( ^3|) are in agreement with the general expression 
for the two lowest-order HTSE expansion terms for C (T) 
of the 5 = 1/2 alternating-exchange Heisenberg chain in 
Eq. ( p^ ) with alternation parameter a = 1. p 

In a later section, the Bethe ansatz C(T) datael will 
be fitted to obtain a function accurately representing the 
C{T) of the 5 = 1/2 AF uniform Heisenberg chain. In 
order that we are not required to change our fitting equa- 
tions from those we use for fitting magnetic susceptibility 
data, the coefficients for the series inverted from that in 



Eq. (43a) are required. We obtain 
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16*2 
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16' 
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(44a) 
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(44b) 




-(C/t)-2/3 

-(C/t)-(2/3)[1+(3/8)g^: 
-(C/t)-(2/3)[1 +(3/8)9^3,9"] 

-(C/t)-(2/3)[1+(3/8)9^a,9'+a,9l 

= 1 .5374(3) 
= 3.125(11 
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FIG. 11. (a) Difference AC/t between the electronic spe- 
cific heat coefficient from the Bethe ansatz data (Ref. ^ and 
the nominal coefficient of 2/3 (top data set), plotted vs re- 
duced temperature t. Moving down the figure, successive data 
sets show the influence of correcting for cumulative logarith- 
mic correction terms, (b) Expanded plots at low tempera- 
tures. The reduced temperature is t — ksT/J and we have 
set N = kB = 1- In both (a) and (b), the lines connecting the 
data points are guides to the eye. 



2. Low-temperature logarithmic corrections 

At first sight, from Fig. || there appear to be no singu- 
larities in the temperature dependence of the specific heat 
for the 5 = 1/2 AF uniform Heisenberg chain. However, 
if the electronic specific heat coefficient C{T)/T is exam- 
ined in detail, one sees anomalous behavior at low tem- 
peratures. Shown as the top curve in Fig. 0(a) is a plot 
of the difference between the electronic specific heat co- 
efficient and its zero-temperature value, /S.C {t) / N k-Qt = 
[C{t) - (2/3)*]/(7VfcB*) for < i < 0.1 [compare with 
Fig. 11(b)] . From this figure, there is still nothing par- 
ticularly strange about the data. However, upon further 
expanding the plot to study the range < i < 0.005 as 
shown in Fig. |l^(b), we see that AC/Nkst is developing 
an infinite slope as i ^ 0. This is the signature of the 
existence of logarithmic corrections to the specific heat 
at temperatures t <C 1, just as it was for the magnetic 
susceptibility. 
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FIG. 12. Log-log plot vs temperature T of the dif- 
ference between our approximate logarithmic expansion 
C'it,iog/(A'^fcB) of Lukyanov's theory (Ref. ^ and his exact 
prediction Cit,g/{NkB) for the low-temperature limit of the 
magnetic specific heat of the spin 5 = 1/2 antiferromag- 
netic uniform Ifeisenberg chain, for the temperature range 
10"^° < t < 0.5. 



Kliimper^ Lukyanov,! and others have found a log- 
arithmic correction to the \ow-t Umit in Eq. (38) 
Lukyanov's exact asymptotic expansion for the free en- 
ergy per spin in zero magnetic field is 



/ 



-J In 2 



3J 



33/2(fcBr)^ 

107rJ3 



(45) 



wher e q(t /to) and to are the same as given in Eqs. ( p5q ) 
and (|35d[ ), respectively, and where g(t) was plotted in 
Fig. gTrhe specific heat at constant volume is calculated 
using C = -Td'^f/dT^, yielding 



Cit.g(r) ^ 2kBT 

NkB 



3J 



1 3 3 
8^ 



This formula shows that the electronic specific heat co- 
efficient C{T)/T increases quadradically with T at low 
T (after subtracting the logarithmic corrections). The 
numerical prefactor of the term is 1.98478 • • •. If the 
approximate expansion for g{C) in Eq. ( ^6|) is substituted 
into Eq. (Hq), one obtains 



Clt,log(J') 

Nkn 



2kBT 
3J 



(2/:)= 



o 



2(3^/2) /fcBT\3 



(47) 



where the prefactor 3/8 in the logarithmic correction 
term was found independently by Kliimpcr,Ll confirming 
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FIG. 13. Semilog plot vs temperature T of the difference 
AC = C — Cit,g (•) between the Bethe ansatz numerical spe- 
cific heat C{T) data (Ref. |^ and Lukyanov's theory (Ref. ^) 
[C'it,g(r)] for the spin 5* = 1/2 antiferromagnetic uniform 
Heisenberg chain at low temperatures. The error bar on each 
data point is an estimated uncertainty in the theory due to 
higher order correction terms that were not calculated. Also 
shown is the deviation AC = C - (2/3)iVfc|r/J (o) of the 
numerical data from the T ^ hmit of C{T). 



Refs. 61 and ^ The difference between Cit.iog(r) and 
C'it,g(7Tis plotted vs temperature in Fig. n2l where the 
difference becomes > 10^^" only for t ^ 10 . 

Shown in Fig. O is the deviation AC /NkB (•) of the 
Bethe ansatz dataO from Lukyanov's theoretical predic- 
tion in Eq. (^). For temperatures t < 1Q~^, the agree- 
ment is better than 10~*. At higher temperatures, the 
uncertainty in the theoretical prediction due to the un- 
known 0{g^) and higher order correction terms becomes 
an important factor in the comparison. The length of 
the error bar on each data point in Fig. ^ has arbitrar- 
ily been set to {4:/3)tg^{t) [d. Eq. (||)]; the 0{g) uncer- 
tainty in the term is negligible compared to this. Also 
plotted in Fig. ^ is the deviation of the numerical data 
from the extrapolated linear low-T behavior (o). A com- 
parison of the two data sets indicates that the 0{g^) log- 
arithmic correction term is responsible for at least most 
of this latter difference for temperatures t ^ 0.001. 

A more rigorous evaluation of the influence of the 
above logarithmic correction term is obtained by correct- 
ing for it in the plot of AC/t vs t, as shown by the sec- 
ond curve from the top in each of Figs. |ll|(a) and |ll|(b). 
From the latter figure, we infer that although subtracting 
this correction term from the data helps to remove the 
zero-temperature singularity, a singularity is still present 
but with reduced amplitude. This means that additional 
logarithmic correction terms are important, within the 
accuracy and precision of the data. Another indication 
of this is shown in Fig. |lj, where we have plotted AC/t^ 
vs t. According to Eq. (p6[), after accounting for the 
logarithmic correction term(s), the result should be in- 
dependent of t at low t. Instead, both before and after 
accounting for the log correction term, there is a strong 
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FIG. 14. Coefficient AC/t^ of the expected dependence 
of the electronic specific heat coefficient at low temperatures 
(top data set). Moving down the figure, successive data sets 
show the influence of correcting for cumulative logarithmic 
correction terms. If all logarithmic corrections were accounted 
for, the data would become independent of t for f — > 0. 
Here, the reduced temperature is i = k^T / J and we have 
set N = = 1. The lines connecting the data points are 
guides to the eye. 

upturn at low temperatures although the strength of the 
upturn is smaller after subtracting the influence of the 
log correction term. I-I 

The numerical Bethe ansatz specific heat datao are 
sufficiently accurate and precise that we can estimate 
the coefficients of the next two logarithmic correction 
(ff''; 9^) terms in Eq. ( ^6| ) from these data as fol- 
lows. From Eq. (Hq), if wc plot the numerical data as 
[C{t)/{NkBt) - (2/3) (1 + 3^3/8)]// vs g at low tem- 
peratures, where the term can be neglected, and fit 
the lowest t data by a straight line, the y intercept for 
5 — > gives the coefficient of the g"^ term and the slope 
gives the coefficient of the term. This plot is given 
in Fig. |l^. This type of plot places extreme demands on 
the accuracy of the data. Even so, we see that the data 
follow the required linear behavior even at the lowest 
temperatures. We fitted a straight line to the data from 
i = 5 X 10~^^ up to a maximum temperature i™'^'^. The 
fit parameters and rms deviation held nearly constant for 
^max ^ 5 X 10-15 (11 data points) up to t^'''' 5 x 10"^ 
(18 data points), but both quantities changed rapidly 
upon further increasing t"^^^. The fit for t™^^ = 5 x 10~^ 
is shown as the straight line in Fig. 



15, along with the 



fit parameters. From the parameters of the fit [after ac- 
counting for the prefactor of 2/3 in Eq. (^)], we include 
our estimated coefficients in Eq. (Efl) explicitly as 



Cit,g(r) ^ 2fcBr 

NkB 3J 
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Fit: 1.0249(2) + 2.083(7) g 
Fit range: 0.017<g< 0.052 
(5x 10""<t<5x 10"^) 
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FIG. 15. Plot showing the estimation of the coefficients of 
the C((/'') and 0{g^) logarithmic correction terms in Eq. ( ^ ) 
for the magnetic specific heat of the 5=1/2 AF uniform 
Heisenberg chain. The reduced temperature is t = k^T / J 
and we have set N = k^ = 1. 



a4 = 1.5374(3) , ag = 3.125(11) 



(48b) 



The infiuences of these g'^ and g^ logarithmic correc- 
tion terms on the data in Figs. |ll| and ^ are shown as 
the two additional data sets in each figure, where ac- 
counting for these two terms is seen to largely remove 
the remaining singular behavior as i — > 0. From Fig. 
we can now estimate that the coefficient of the term 
in Eq. ( ^ ) is a little larger than 2, contrary to the ex- 
act value 1.98478- • •. The magnitude of this difference is 
about as expected from the 0{g) logarithmic correction 
to the term, since g(t ~ 0.1) ~ 0.1. The remaining 
upturn at low temperatures in Fig. ^ is due to residual 
logarithmic corrections which are not accounted for. 

If the approxi mate expansion for g{C) in Eq. ( |3^ ) is 
inserted into Eq. (48a), one obtains 



2fcBT 



1 



3J {2Cf 
ln/:[181n(£) - Ua^ 



91n(£) - 1604 
(2Z)4 
9] + 32a5 



{2Cf 



O 



(2£)6 
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IV. FITS TO x'(t) AND C{t) OF HEISENBERG 
SPIN LATTICES 

A. General X*(*) At considerations 

The general expression we use to fit theoretical numer- 
ical x*(t) data for 5* = 1/2 Heisenberg spin lattices is 



e-^tit/* / ^ 



(50a) 
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(50b) 



(<?) 



where the orders q and r of the Pade approximant P^^,^ 
are often constrained by the behavior of x*{t) low t, 
and the fitted gap Ag^ is not necessarily the same as 
the true gap. At high t, x*{t) in Eqs. ( |50| ) approaches 
the Curie law l/(4i) as required [for a gen eral spin 
lattice, the numerical prefactor 1/4 in Eq. (50a) would 
be replaced by S{S +l)/3]. 

The Nn and Z)„ parameters in Eq. ( |50b| ) are not in 
general independent if one or more of the HTSE condi- 
tions in Eqs. (5b) and (j5^ are invoked. For example, for 
n = 1-3 one finds 



Di = {di+Ni)-A 



fit ) 



D2 = (da + diNi + N2) ~ A*M + Ni) + 



^fit 



(51a) 



(51b) 



^3 = (ds + d2Ni + diN2 + N3) - A*M + diNi + N2) 



-^{di 



In general, one has 

n n—p 
p—O m—O 



Ni) 



(51c) 



(-Al)P 

pi 



(52) 



A fit of experimental or theoretical x*(0 data by 
Eqs. (pG) can be constrained b y in serting one or more 
of Eqs. (|51|) and ( p2|) into Eq. (50b). These constraints 
are especially useful for high-t extrapolations when x* (t) 
data are not available for high temperatures t ^ 1, 
and/or to reduce the number of fitting parameters re- 
quired to obtain a fit of specified precision. In the fol- 
lowing fits to the numerical x*(^) data for the dimer, the 
uniform chain, and finally our QMC and TMRG data for 
the alternating-exchange chain, the three constraints in 
Eqs. ( |5l| ) on Di, 1)2, and D3, respectively, are enforced 
in each case, where di, c?2, and ^3 for the alternating- 
exchange chain are given in Eq. (p^). 

All of the fits reported in this paper were carried 
out on a 400 MHz Macintosh G3 (B&W) computer 
with 1GB of RAM. Most fits were implemented us- 
ing the program Mathematica 3.0, although a few of 
the simpler ones (fits to experimental data) were done 
using KaleidaGraph 3.08c. The fits using Mathematica 
sometimes required prodigious amounts of memory, e.g., 
930 MB for the 28-parameter fit to the combined 2551 
data point QMC and TMRG x*ia,t) data set for the 
alternating-exchange chain in Sec. ^ below. 



B. Fit to X* {t) of the S = 1/2 antiferromagnetic 
Heisenberg dimer 

The spin gap of the 5 = 1/2 Heisenberg dimer is A = 
J, where J is the antiferromagnetic exchange constant 
within the dimer. The spin susceptibility and its low- 
temperature limit are given by Eqs. (|^). The x*(0 
plotted in Fig. 0(a) for 0.02 < t < 4.99. In order to 
later obtain a continuous fit function for x*i'^i ^) fo'' the 
entire range < a < 1 of the alternating-exchange chain, 
it is necessary to first obtain a high accuracy fit to the 
exact expression (p^ for the dimer by our general fitting 
function in Eqs. (pG), in addition to Fit 1 obtained for 
the uniform chain below. The form of our fit function in 
Eqs. ( po|) allows both the low- and high-t limiting forms of 
X*{t) for the dimer to be exactly reproduced. The low-t 
limit in Eq. ( ^ ) requires that r — q and that Dq — iVg/4 

in the Fade approximant T^l^j' take Agt = A, so 

the total number of fitting parameters is 2q — 4. 

We fitted the 498-point double-precision representa- 
tion of x*{t) in Fig. |l|(a) from t = 0.02 to t = 4.99 by 
Eqs. (|5^) using the above constraints. The variances of 
the four fits for q = r = 4, 5, 6 and 7 were 2.5 x 10""'^^, 
1.17 X 10~l^ 5.3 X 10"", and 5.6 x lO'^^ respectively, 
showing that Eqs. ( |50| ) have the potential for very high 
accuracy fits with a relatively small number of fitting pa- 
rameters. The six Nn {n = 1-5) and parameters of 
the fit for q, r = 5 are given in Table |, along with Di, D2, 
and D3 computed from Eqs. (^ij) and D5 — Nc,/4,. The 
Fade approximant T'lgj in the fit function has no poles 
or zeros on the positive t axis. The fit is shown by the 
solid curve in Fig. |l|(a), and the deviation of the fit from 
the exact susceptibility in Eq. ( pa|) is plotted versus t in 
Fig. %{h). 



C. Fits to X*(*) of the S = 1/2 Antiferromagnetic 
Uniform Heisenberg Chain 

Fit 1: 0.01 < t < b. Fits to the uniform chain 
X*{i) calculated by Eggert, Affieck and Takahashi3 fca: 
limited temperature regions were obtained previously.u 
Here we obtain a fit (Fit LI to the higher accuracy data 
of Kliimper and JohnstonO for the temperature region 
0.01 < i < 5 (999 data points) using Eqs. (|^), the re- 
sults of which will be utilized later in the fit function for 
t > 0.01 for our QMC and TMRG alternating-exchange 
chain x*(ct,^) data. This uniform chain fit can be accu- 
rately extrapolated to arbitrarily high t. 

The requirement that x*(^ — > 0) is a finite non-zero 
value requires A^^ = and r = g -I- 1 in Eqs. (p0|). We 
found that using q = 5 and r = 6 produces a fit suffi- 
ciently accurate for use in the fit function for our QMC 
and TMRG calculations for the alternating chain. The 
seven iV„ {n = 1-5) and £)„ {n = 4-6) parameters ob- 
tained for the fit with q = 5, r = 6 are given in the 
column labeled "Fit 1" in Table |, along with Z?i, D2, 
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The quality of Fit 1 does not approach the limitation 
imposed by the absolute accuracy of the data (1 x 10^^). 
For an ideal fit, the variance is expected to be ~ 10~^^ 
and the relative rms deviation ~ 0.01 ppm. As can be 
inferred from Fit 2 in the following section, the reason 
that Fit 1 cannot be optimized to this extent is due to 
the t — critical point and associated logarithmic di- 
vergence in the slope of x*(0 as t ^ 0; this divergence 
cannot be fitted accurately by a finite polynomial or Fade 
approximant. We attempted to improve the accuracy of 
the fit over the same temperature range 0.01 < t < 5 by 
replacing the Fade approximant 'Pj^^j in the fit function 

by 'Pjyj , which incorporates two additional fitting param- 
eters. The variance improved somewhat to 2.18 x lO"""^^ 
and the relative rms deviation improved slightly to 12.2 
ppm, but the Fade approximant developed a pole at 
1/t = 129.23, and hence this fit was discarded. Although 
the temperature at which this pole occurs is below the 
fitted temperature range, as a general rule we cannot al- 
low poles in the fit function at low temperatures because 
of problems that can occur when using the fit function 
to model experimental data which include data at tem- 
peratures lying below the fitted temperature range of the 
fit function. In fact, we will encounter this situation fre- 
quently in modeling experimental data later. For exam- 

TABLE I. Fitted parameters for x* (^) of the 5* = 1/2 antiferromagnetically coupled Heisenberg dimer (a = 0) [Eqs. ( |5c| ) 
with = 1] and x*(i) and Cit) [Eqs. @] for the uniform chain (a = 1). X*(*) Fit 1 for the uniform chain (0.01 < t < 5) 
[Eqs. (^) with AJt — 0] uses powers of 1/t only, whereas x*(*) Fit 2 (0 < t < 5) [Eqs. (^)] also incorporates logarithmic 
correction terms. 



parameter 


X*(a = 0) 


X*(a = 1) Fit 1 


X*{a = 1) Fit 2 


C{a = 1) 


iVi 


0.6342798982 


-0.053837836 


-0.240262331211 


-0.018890951426 




0.1877696166 


0.097401365 


0.451187371598 


0.024709724025 


iVs 


0.03360361730 


0.014467437 


0.0125889356883 


-0.0037086264240 




0.003861106893 


0.0013925193 


0.0357903808997 


0.0030159975962 




0.0002733142974 


0.00011393434 


0.00801840643283 


-0.00037946929995 


iVe 






0.00182319434072 


0.000042683893992 


N7 






0.0000533189078137 




Ns 






0.000184003448334 










1.423476309767 










0.341607132329 




ti 






5.696020642244 






-0.1157201018 


0.44616216 


0.259737668789 


-0.51889095143 


D2 


0.08705969295 


0.32048245 


0.581056205993 


0.59657583453 


Da 


0.005631366688 


0.13304199 


0.261450372018 


-0.15117343936 


Di 


0.001040886574 


0.037184126 


0.142680453011 


0.074445241148 


D5 


0.00006832857434 


0.0028136088 


0.0572246926066 


-0.0024804135233 


De 




0.00026467628 


0.0176410851952 


-0.00053631174698 


D7 






0.00390435823809 


0.00082005310111 


Ds 






0.000119767935536 


-0.00010820401214 


I?9 








0.000011991365422 


ai 








-0.000015933393 
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FIG. 16. Semilog plot vs temperature T of the deviation of 
Fit 1 from the magnetic susceptibility calculations of Kliimper 
and Johnston (Ref. ^ for the 5* = 1/2 uniform antiferromag- 
netic Heisenberg chain. 

and Z?3 computed from Eqs. (pl|). The Fade approxi- 
mant P^gj^ in the fit function has no poles or zeros on the 
positive t axis. The deviation of the fit from the data is 
plotted in Fig. |l^ The variance of the fit is 2.97 x 10"^^ 
and the relative rms deviation of the fit from the data in 
the fitted t region is 14.5 ppm. Extrapolation of the fit 
to higher temperatures is very accurate. 
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pie, for NaV205, t = 0.01 corresponds to an absolute 
temperature T « 7K, whereas the experimental data 
and modeling extend down to « 2 K. 

Fit 2: < t < 5. We can greatly improve the accu- 
racy of the fit compared to that of Fit 1, and extend the 
fit to t = 0, by restricting the high-temperature limit of 
the fit and using in the fit function one or more low- 
temperature logarithmic correction terms discussed in 
Sec. [II A 2. In particular, in this section we obtain a 
very high precision fit (Fit 2) to the exact t = O^value 
and to the calculations of Kliimper and JohnstonQ over 
the entire temperature range 5 x 10~^^ < t < 5 of the cal- 
culations. We do not use this fit in our formulation of the 
fit function for the alternating-exchange chain. However, 
Fit 2 will be generally useful for evaluating the accuracy 
of other theoretical calculations of x*{t) for the uniform 
chain, such as our TMRG calculations to be presented 
below, and for modeling appropriate experimental x{T) 
data whose scaled upper temperature limit is below t — 5. 

We initially formulated a fit function utilizing a mod- 
ified Fade approximant in which the last term of the 
numerator and/or denominator contained the Xit log 
pansion in Eq. (|3^), such that the low-temperature ex- 
pansion of the fit function yielded Xit log ^o lowest or- 
ders in t. The best fit to the data from t = 5 x 10 
to 2.5 (777 data points) was unsatisfactory, with a vari- 
ance V ~ 2.4 X 10~^^ and a relative rms deviation 
Crms = 45 ppm. Allowing an arbitrary coefficient in 
place of the exact value a/S/tt^ yielded an improved fit 
with V — 1.1 X 10~^^ and CTrms = 9.6 ppm. However, 
this fit was still unsatisfactory, given the high absolute 
accuracy (1 x 10~^) of the data. From these results it 
became clear that a fit function which can fit the data 
to much higher accuracy over such a large temperature 
range would indeed have to include an expression Xi*og(i) 
containing logarithmic correction terms, but where the 
form and/or coefficient of one or more of these terms 
would have to be empirically determined by trial and 
error. This process yielded the formulation we now de- 
scribe. 

The X*og(^) function is incorporated into our fit func- 
tion in Eqs. ( [50| ) as follows. As in Fit 1, the finite value 
of x*(0) requires A^t 
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in Eq. ( pOaj ) and 



in the Fade approximant T'l^' (t) in Eq. ( ^0b| ) 

two terms highest order in 1/t in 'P\y^ (t) (one each in the 
numerator and denominator) dominate the fit as t —> 
and become small for t ^ 1, relative to the other terms 
in the numerator and denominator, respectively, we in- 
corporate X*og(^) info the last term in the numerator of 

a modified 'P'}y^{t). Trial fits showed that to obtain the 
optimum accuracy of the fit required q — 8 and r = 9. 
Our final fit function for Fit 2 is 

i + [ELi^«Al+4iV8xrog(0A' 
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FIG. 17. (a) Semilog plot of the absolute deviation of Fit 2 
(0 < t < 4) from the magnetic susceptibility calculations of 
Kliimper and Johnston (Ref. ^) versus temperature T. (b) 
Expanded plot of the data in (a) at the higher temperatures. 
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C = \n{h/t) , 



(53b) 



(53c) 



(53a) 



subject to the three constraints on Z?i, D2, and Da in 
Eqs. ( ^ which are required by the HT SE. Two of the 
four logarithmic correction terms in E q. ( 53b| ) are identi- 
cal to the first two such terms in Eq. (34a). By construc- 
tion, the exact x*(0) = I/tt^ is fitted exactlv. 

We fitted all of the numerical x*(^) datajj calculated 
over the range 5 x 10^^^ < ^ < 5 (1119 data points), 
by Eqs. (^3|). The 19 fitting parameters of the fit func- 
tion (53|), which are Nn in = 1-8), Dn {n = 4-8), 
-^81, ^^82 and ti, are given in the column labeled "Fit 2" 
in Table ^, along with Di, D2, and computed from 
Eqs. (^l|). The data to parameter ratio is 59. Th e den om- 
inator of the modified Fade approximant in Eq. ( |53a| ) has 
no zeros for any real positive t. The fit is shown in the 
low-temperature region < t < 0.02 in Fig. |(b) [ over 
the larger t range plotted in Fig. ^(a), the fit is indis- 
tinguishable from the data and is therefore not plotted 
there] . 
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The deviation of Fit 2 from the numerical data for 
10"^^ < t < 5 is plotted vs log^oi in Fig. ^(a), and 
an expanded plot at the higher temperatures is shown in 
Fig. |l3(b). Due to a logarithmic divergence in Xiogi'^) 
t = ti = 5.696, Fit 2 should not be used (e.g., for mod- 
eling experimental data) at temperatures i 5. The 
variance of the fit is 9.8 x 10~^^, and the relative rms de- 
viation is cTrnis — 0.087 ppm. These values are both much 
smaller than for Fit 1 above. The relatively large number 
of fitting parameters in Fit 2 is justified a posteriori by 
the extremely high quality of the fit over a temperature 
range spanning 25 orders of magnitude. 



D. Fit to C{t) for the S = 1/2 antiferromagnetic 
uniform Heisenberg chain 



The logarithmic corrections to the magnetic spe- 
cific heat C {t) at low temperatures, discussed above in 
Sec. [II B 2, do not pose as serious a problem for fitting 
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the data as for x*(0; because the strength of these log 
corrections is much smaller for C{t) than for 'X*{t). In 
addition, since here we fit C(t), and not the electronic 
specific heat coefficient C{t)/t, the influence of the log 
corrections is ameliorated by the multiplicative leading 
order dependence of C{t). Even so, in order to obtain 
the otjtimum fit to the highly accurate Bethe ansatz C (t) 
datajj we found it necessary to take the influence of the 
logarithmic corrections into account. j-. 

Our flt to the Bethe ansatz C{t) dataJj some of which 
were shown previously in Fig. @(a), was carried out in 
two stages. First, the data from t = 0.01 to the maxi- 
mum temperature t — 5 of the calculations were fltted 
by the Fade approximant vjy^ in Eq. ( ^Ob| ) with a pref- 
actor 3/(16^2) to satisfy the HTSE in Eqs. (||) to lowest 
order in 1/t. The orders q and r of T'^'-.^j were chosen to 
satisfy r — q + 3 so that C{t — > 0) cx i. To obtain a 
fit of the required accuracy (see the fit deviations given 
below) we found that g = 6 and r = 9 are of sufficiently 
high order. Due to the presence of the log corrections at 
very low i, we did not require the parameters Nq and Dg 
to yield the exact coefficient 7 = 2/3 in the expression 
C{t)/Nk-B — 7^, in a \ow-t expansion of the fit function. 
We also found that to obtain the best fit, o nly t he one 
additional HTSE constraint (on Di) in Eq. (51a) (with 
A^j = 0) could be used. It was quite difficult to find 
the region in parameter space in which the absolute min- 
imum in the variance of the fit resided; the initial starting 
parameters usually flowed to regions with local variance 
minima in them with much larger values (by two to four 
orders of magnitude) than the smallest variance we 
ultimately found. Then the deviation of the fit from all 
the data for 5 x 10^^^ < t < 5 was computed. The fit de- 
viations for t > 0.01 were very small [0(10"^)], but the 
log corrections which become most important at lower 
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FIG. 18. Semilog plot of the absolute deviation of the fit 
(10~^^ < ksT / J < 5) from the specific heat calculations of 
Kliimper and Johnston (Ref. ^ in Fig. ^(a) versus temper- 
ature T. The fit deviation is negligible at low temperatures 
T < 10~®J/feB- The lines connecting the data points are 
guides to the eye. 

temperatures resulted in fit deviations aX t < 0.01 an 
order of magnitude larger. We therefore fitted the fit de- 
viation versus t for 0<i<0.1bya separate empirically 
determined function Fit), so the net fit function consists 
of the Fade approximant fit function minus the fit func- 
tion to the low-t fit deviations. In the final fitting cycles 
the two functions were refined simultaneously. 

Our final fit function for C{t) in the range < i < 5 is 



C{t) 



V\fAt)-F{t) 



NkB 16t2 (9) 



F{t) — sin 



27r 



a2 + cist 



05^ e' 



-net 



(54a) 



(54b) 



(54c) 



subject to the constraint on Di in Eq. (51a) which 
is required by the HTSE. By construction, the exact 
C(0) = is fitted exactly. The 20 fitting parameters, 
Nn {n = 1 — 6), Dn (n = 2 — 9) and a„ (n = 1 — 6), are 
given in Table |, together with the constrained parameter 
Di computed from Eq. (51a) with A^^. = and di given 
in Eq. (44b). The deviation of the fit from the data is 
shown in a semilog plot vs temperature in Fig. |l^. The 
maximum deviations of ~ ±4 x 10~* occur at t « 0.3. 
The absolute rms deviation of the fit from all the data 
(1119 data points), which extend over the temperature 
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range 5 x lO'^s < t < 5, is 1.34 x 10"®, and the rela- 
tive rms deviation for 0.01 < t < 5 (999 data points) is 
0.50 ppm. 

At high temperatures, our C{t) fit function reduces 
by construction to the lowest order 1/P and terms 
of the HTSE of C{t) in Eqs. (||), so extrapolation of 
our C{t) fit function to arbitrarily higher temperatures 
should be very accurate (see Fig. 18). In particular, even 
though our fit was to C(t) and hence not optimized as 
a fit to the electronic specific heat coefficient C{t)/t, the 
magnetic entropy S at t — oo computed from our C(t) 
fit function is 



S{t 



Nkj: 



C{i) 
Nknt 



dt = 0.693147 235 



(55) 



which is the same as the exact value ln2 = 0.693 147 181 
to within 8 parts in 10*. This agreement reflects well on 
our fit function, and of course also strongly canfirms the 
high accuracy of the Bethe ansatz C(t) data.l3 



E. Fit function for the S = 1/2 AF 
alternating-exchange Heisenberg chain 

Here we formulate a single two-dimensional (a, t) func- 
tion to accurately fit numerical calculations of %* {a, t) for 
the 5=1/2 alternating-exchange Heisenberg chain for 
the entire range < a < 1, and for the entire tempera- 
ture range t > 0.01 over which our Fit 1 for x*(0 of the 
uniform chain is most accurate, subject to four general 
requirements as follows, (i) The HTSE of the x*(q^i^) 
fit function must give the correct result to 0(1/^^), as 
satisfied by the fit functions for the dimer and uniform 
chain (Fit 1) susceptibilities above, so that the fit can be 
accurately extrapolated to higher temperatures, (n) We 
require the x*(a,t) fit function to become identical with 
those found above for the isolated dimer and for the uni- 
form chain (Fit 1) when a — and a = 1, respectively. 



As discussed above in Sec. [I C, at any finite temperature, 
X*{S,t) in the variables S and t is an even (analytic) func- 
tion of 6. Therefore, as a minimum accommodation of 
this fact, (iii) we require that the fit function for x*('^i t), 
when transformed to the form x*{^,'t), must have the 
property dx*{S,t)/dS\s=Q = at all finite temperatures. 
This requirement is clearly the minimum necessary in or- 
der to accurately interpolate the fit vs a for a — s- 1 at 
each t, and to thereby accurately model the susceptibility 
of materials which are in or near this limit. Finally, the 
QMC and TMRG calculations of x*(q;, t) to be presented 
below are sufficiently accurate and cover sufficiently large 
ranges of a and t with sufficient resolution that (iv) we 
require the nonanalytic energy gap A(a) [see Eqs. ( |l9| ) 
and (po|)] to be included in the fit function in order to fit 
the data for a < 1 at i ^ 1, so as to avoid the alternate 
necessity of including high-order power series in a and t 
in the fit function. We note that according to Eq. ( [I9b| ) 
or (|2^), dA*{S)/dS\s=o = oo. The major obstacle we 



faced in formulating the fit function for x*(a,f) was to 
simultaneously satisfy both requirements (iii) and (iv), 
which at first sight seem to require mutually exclusive 
forms for the fit function. 

We found that these four requirements can all be sat- 
isfied by an extension of the form of the fit function in 
Eqs. (|g) which was used above for the isolated dimer 
and for the uniform chain Fit 1. This extension consists 
of using a modified Fade approximant 'Pm|g| in the fit 



function in place of the former vj"^} . The fit function is 



X*{o-,t) 



e 



-Af*t(a)/t 



At 



(r)- 



(7) 



(a,t) 



(56a) 



ELo ^»/t"] + [Nna + N72a^)iA*Jt)y /t^ 
ELo Du/t-] + {Dsia + D82a2)(A5/i)^e(^o*-AS.)/7i8 ' 

(56b) 



A*,(a) = l-ia-2a2 + ^a3 ^ 



A*(a) = (l-a)3/4(l + a)i/4 

-|- 5iq(1 — a) + g20?{^ — otf' 

4 

iV„(a) = ^ 7V„„a™ (n = l-6), 

m=0 



Dnia) 



(n = 1 - 7) 



(56c) 



(56d) 



(56e) 



(56f) 



(56g) 



To satisfy requirement (i), Di{a), D2{a), and D^i^a) 
are determined from the Ni{a), N2{a), and N^^a) fitting 
parameters according to the three constraints in Eqs. ( |5l| ) 
demanded by the HTSE. In order to satisfy requirement 
(ii), the {Nno, Dno} parameters are set to be identical 
with those determined above for the dimer, and we re- 
quire {Nn{l), D„(l)} to be identical with the correspond- 
ing fit parameters determined above in Fit 1 for the uni- 
form chain. In order to satisfy requirement (iii), the Nnm 
and Dnm coefficients must satisfy 



(n - 2m) {Nnm or £)„„) = , 



(57) 



so that no term appears in the Taylor series expan- 
sions in 6 of the transformed {7V„((5), Dn{S)}. These var- 
ious constraints on the {Nmn, Dnm} parameters reduce 
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the number of independent fitting parameters within 
this set from 50 to 20. Together with the parameters 
^^^71, JY22, A^si, ^82, z in Eq. (p6b|) and 51, 52 in 



Eq. ( 56d ), the total nmnbcr of independent fitting pa- 
rameters in the fit function is 28. 

The quantity A^^(q!) in the exponential prefactor to 

'Pm(gj in Eq. ( ^Gaj ) cannot be set equal to the true non- 
analytic gap A* (a), because this prefactor affects the fit 
at all i, and would not allow requirement (iii) above to be 
fulfilled. In addition, the nonanalytic critical behavior of 
A*(a — > 1) in practice only becomes manifest in x*(ci^i 
at low temperatures t Therefore, we separated the 

spin gap into an analytic part Agj(a) which go es in to the 
argument of the exponential prefactor in Eq. (56a), and 
a nonanalytic part AQ(a) [satisfying requirement (iv)] 
which is placed into the argument of the exponential in 
the last term of the denominator of 'Pin(gj in Eq. (561:) 
and which therefore only becomes important at low tem- 
peratures. The first two terms of Ag^(a) (to order a^) 



in Eq. (56c) are the first two terms of the exact dimer 
series expansion up to O(a^) given by Barnes, Riera and 
TennantE3 for the AF alternating-exchange chain, and 
the last two are included so that d^\^{S)ld&\^=Q ~ 0, 
in accordance with requirement (iii). The nonan alytic 
Ag(a) in Eq. (56d) contains the behavi(p,in Eq. (19a) 
proposed by Barnes, Riera, and Tennant,t2l plus two an- 
alytic terms which are included to adjust the a depen- 
dence for a — !■ 1 but which make no contribution at a = 
or a = 1. Provided that the inequali ty y, 2; > 4/3 is sat- 
isfied by the powers y and z in Eq. (56b), the last term 



in each of the numerator and denominator of T'm^gj (a, t), 

when transformed to the variables (^, i), has a partial 
derivative with respect to b which is zero at (5 = 0. 

We have now shown that at (5 = (a = 1), the par- 
tial derivative of each part of x* (^, *) with respect to (5 
is zero (if y, z > 4/3, which is confirmed in the actual 
fit later). Hence, the entire fit function has the prop- 
erty 9x*((J, i)/9(5|5=o = at all finite temperatures, thus 
satisfying requirement (iii), despite the fact that the fit 
function contains the nonanalytic l^^ip) as required by 
requirement (iv). 

At the lowest temperatures, the last term in each of 
the numerator and denominator of 'Pm(g) in Eq. (56t) 
should dominate the fit, together with the exponential 

(7) \ V 

prefactor to 'Pm(g) in Eq. (|56a|), so in this limit our fit 
function for < a < 1 becomes 



X*(a,t^O) 



iVyiQ! -I- N'j20? 



A* (a) 



-A'o{a)/t 



(58) 



This expression has the form of Eq. ( |lOa| ) (with ^ = y—z) 
as required in the \ow-t limit. In fact, the forms of the 
last term in each of the numerator and denominator of 
VraiJ) were designed to result in the form of Eq. (10a) in 



the low-t limit, with Aq and t entering the prefactor only 



as their ratio as in Eq. (p^), in addition to being consis- 
tent with requirements (iii) and (iv). One might expect 
the fitted y and z powers to satisfy y — z — ^ = 1 /2 as in 
Eq. (|l0|). However, if a fit of x*(i) data by Eq. (|l0|) is 
not carried out completely within the low-i limit, an ef- 
fective exponent 7 ~ 1 is often inferred [see, e.g., Eq. ([70| ) 
and subsequent discussion, and Fig. ^ below] . Similarly, 
since many of our calculated x*(q;, t) data sets for differ- 
ent a in the fitted temperature range t > 0.01 are not, 
or do not contain extensive data, in the low-t limit, we 
did not impose the constraint y — z = 1/2. On the basis 
of the above discussion we expect the actual fitted values 
of y and z to yield ?/ — z ~ 1. In fact, as will be seen 
in the next section, our fitted parameters y and z give 
y - z = 1.14. 



V. QMC AND TMRG x*(«,T) CALCULATIONS 
AND FIT FOR THE S = 1/2 AF 
ALTERNATING-EXCHANGE HEISENBERG 
CHAIN 

QMC simulations of x* («, t) were carried out on S* = 
1/2 alternating-exchange chains containing 100 spins for 
a = 0.05, 0.1, 0.15, 0.9, 0.92, 0.94, 0.96, 0.97, 
0.98, and 0.99 in various temperature ranges spanned by 
0.01 < t < 4. 

Complementary TMRG calculations of x*{'^^ t) oi S = 
1/2 alternating-exchange chains were carried out for a = 
0.80, 0.82, 0.96, 0.97, 0.98, 0.99, 0.995 and 1, where 
the number of states kept was m = 150 or 256. The 
calculations were carried out for temperatures given by 
l/i = 0.1, 0.2, (1/i)'"'^'^, with (l/i)"'''''< 500 increas- 
ing with increasing a, and comprised a total of 22 370 
(a, t) parameter combinations. The details of the calcu- 
lational method are given in Refs. and |ll|. It should 
be noted that the TMRG calculations by their nature are 
explicitly in the thermodynamic limit. 

The reason for doing TMRG calculations for the uni- 
form chain (a = 1) was to enable comparison of the re- 
sults with the valuesQ computed with the Bethe ansatz 
which have a high absolute accuracy of 1 x 10^^. This 
comparison was done using the above very accurate and 
precise Fit 2 for the Bethe ansatz data. The relative devi- 
ation of the TMRG data from Fit 2 is shown in Fig. |9|(a), 
and an expanded plot for the higher temperature region 
t > 0.01 is shown in Fig. |l|(b). This comparison in- 
dicates that the accuracy of the TMRG calculations for 
both m = 150 and 256 in the range t > 0.01 is-|better 
than 0.1%, which is the same as the estimateE£l made 
previously for m = 80. However, the accuracy of these 
calculations deteriorates rapidly at lower t, to about 3 % 
at the lowest temperatures t « 0.002 for m — 150. 

Since the TMRG calculations extend close to the t = 
limit for most of the above-stated a values, the spin gaps 
can be estimated from these data. Comparisons with 
previous work can then be made of the dependence 
of the spin gap on a. An important question, not yet 
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FIG. 19. (a) Semilog plot vs temperature T of the relative 
deviation of the magnetic susceptibility x of the S = 1/2 
antiferromagnetic uniform spin Heisenberg chain calculated 
with the TMRG technique from that calculated (Ref. |^ using 
the Bethe ansatz. (b) Expanded plot of the data in (a) at the 
higher temperatures. 

answered in previous work, is the approximate a value at 
which the asymptotic critical region is entered upon ap- 
proaching the uniform limit. Performing these estimates 
and comparisons will be postponed to the following sec- 
tions. In the present section, we present the QMC and 
TMRG x*{ct, t) data and obtain a fit to these combined 
data by the fit function formulated in the previous sec- 
tion. 

Some of the results for t < 2 are shown as the filled 
symbols without error bars in Fig. ^0|(a) (the error bars 
are smaller than the data symbols); an expanded plot of 
data for t < 0.4 is shown in Fig. ^(b). [A log-log plot 
of the TMRG x*(Q!,t) data at low t is shown below in 
Fig. Also shown in both figures as the two bounding 
solid curves with no data points are the fits we obtained 
above to x*(i) for the dimer and uniform chain (Fit 1), 
respectively. The data points plotted for a given a value 
are the subset below the upper temperature limits of the 
figures, of the subset of available data points which were 
fitted by our fit function as described below. 

We fitted a combined QMC and TMRG x*(a,<) data 
set containing 2551 selected data points over the temper- 
ature range 0.01 < t < 10. The 802 QMC data points 
covered the ranges 0.01 < t < A and 0.05 < a < 0.99. 



0.20 




FIG. 20. (a) Magnetic susceptibility x versus tempera- 
ture T for the spin S — 1/2 antiferromagnetic (AF) alter- 
nating-exchange Heisenberg chain with alternation parameter 
a = J2/J1 from to 1, as shown. The small filled circles are a 
selection of the calculated QMC and TMRG data, where for 
clarity only a small subset of the available data are plotted. 
The set of curves through the data is obtained from the global 
two-dimensional {a,t) fit function in Eqs. js^ ) with parame- 
ters given in Table |ll The solid curves for a — and 1 are 
plots of the fit function for the dimer and uniform chain sus- 
ceptibilities, respectively, for which no data are plotted. The 
parameter Ji is the larger of the two alternating exchange 
constants, (b) Expanded plot of the fit for a selection of data 
at low temperatures. Error bars are plotted with the QMC 
data in (b), but are not plotted in (a) because they are not 
visible on the scale of this figure. 

The average estimated absolute accuracy of these QMC 
data is 1.7 x 10"''. The best estimated accuracy among 
these QMC data is 7.7 x 10^^ and the worst is 1.5 x 
10"'^, with the better accuracies occurring at the highest 
temperatures. The 1749 TMRG data points covered the 
ranges 0.01 < t < 10 and 0.8 < a < 0.995. We did not 
use all 22 370 TMRG data points in the available data set, 
because this would have weighted the region a ^ 1 too 
heavily in the fit, and in any case a large fraction of these 
are for temperatures below our low-temperature fitting 
limit of t = 0.01. We used the low-temperature data 
to determine the spin gaps as described in the following 
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section. 

We fitted this x*(a,t) data set by Eqs. with the 
constraints on the parameters discussed above. Obtain- 
ing a rehable 28-parameter two-dimensional fit to these 
data over the full above-cited ranges of t and a, with 
no poles in the fit, posed a very difficult challenge. The 
particular choice of starting parameters and the detailed 
sequence of refinements were found to be important to 
avoiding poles in the final fit. Since there are a total of 
28 parameters in the fit function for 2551 data points, 
the data to parameter ratio is 91. The number of fitting 
parameters seems large, until it is realized that we are si- 
multaneously fitting x*it) data for 29 different a values, 
so on average a x*(^) data set for a given a value is fit- 
ted by a single parameter. A weighting function was not 
included during the variance minimization, because we 
were interested in obtaining a fit which treated all the 
data points the same on an absolute scale; this choice 
optimizes the fit for use in modeling experimental data. 

The parameters of the fit are given in Table ||, where 
we have also included the constrained parameters for 
completeness and for ease of implementation of our fit 
function by the reader. From Eqs. (pT|), the constrained 
parameters D2 and D3 contain products of the third- 
order (in a) polynomial Ag^ with itself and/or with 
the fourth-order A^i fitting polynomial, so D2 and D3 
are of seventh and tenth-order, respectively. The two- 
dimensional fit is shown as the set of solid curves in 
Fig. The variance of the fit is w = 3.77 x 10"^ 
The absolute rms deviation ^/v « 1.9 x 10^* is about 
the same as the average estimated accuracy of the QMC 
data noted above, indicating that the fit function is ap- 
propriate and that the fit is a reliable representation of 

TABLE II. Parameters in the fit function [Eqs. (B^] for x*{a,t) of the S = 1/2 antiferromagnetic alternating-exchange 



Heisenberg chain. Note that D2 and D3 


are respectively of seventh and tenth order 


in a. 






parameter 


m = 


m — 1 


m = 2 




m = 3 


m = 4 


iVlm 


0.63427990 


-2.06777217 


-0.70972219 




4.89720885 


-2.80783223 


iV2m 


0.18776962 


-2.84847225 


5.96899688 




-3.85145137 


0.64055849 


Nam 


0.033603617 


-0.757981757 


4.137970390 




-6.100241386 


2.701116573 


N4,n 


0.0038611069 


0.5750352896 


-2.3359243110 




2.934083364 


-1.1756629304 


A5,„ 


0.00027331430 


-0.10724895512 


0.40345647304 




-0.48608843641 


0.18972153852 


N(im 





0.00578123759 


-0.02313572892 


0.02892774508 


-0.01157325374 


Arm 




2.59870347x10"^ 


-2.39236193x10" 


-7 






Dim 


-0.11572010 


-1.31777217 


1.29027781 




3.39720885 


-2.80783223 


D2m 


0.08705969 


-1.44693321 


5.09401919 




-10.51861382 


8.97655318 




5.75312680 (m = 5) 


-11.83647774 (m = 6) 


4.21174835 (m = 


7) 






Dzm 


0.00563137 


0.65986015 


-1.38069533 




-0.09849603 


7.54214913 




-22.31810507 (m = 5) 


27.60773633 (m = 6) 


-6.39966673 (m = 


7) 








-15.69691721 (m = 8) 


13.37035665 (m = 9) 


-3.15881126 (m = 


10) 






Dim 


0.0010408866 


0.1008789796 


-0.9188446197 




1.6052570070 


-0.7511481272 


D^m 


0.0000683286 


-0.1410232710 


0.6939435034 




-0.9608700949 


0.4106951428 


Dfim 





0.0367159872 


-0.1540749976 




0.1982667100 


-0.0806430233 


Djm 





-0.00314381636 


0.01140642324 




-0.01338139741 


0.00511879053 


Dsm 




1.25124679x10"'^ 


-1.03824523x10" 


-7 








92 


y 


z 








0.38658545 


-0.20727806 


4.69918784 


3.55692695 
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FIG. 21. Deviation in absolute units of the fit function for 
the magnetic susceptibility x versus temperature T for the 
spin S = 1/2 antiferromagnetic (AF) alternating-exchange 
Heisenberg chain, with alternation parameter q = J2/J1 
from 0.05 to 0.995, from the QMC data (a) and TMRG data 
(b). The parameter Ji is the larger of the two alternating 
exchange constants. The absolute rms deviations of the re- 
spective data from the fit are given in the figures. The fit 
function is given in Eqs. (p^) with the parameters in Table pi 
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FIG. 22. Deviation in absolute units of the fit function for 
tlie magnetic susceptibility x* versus temperature T for the 
spin 5=1/2 antiferromagnetic (AF) alternating-exchange 
Heisenberg chain, from L — 400 (•) and L = 800 (open 
squares) QMC data for a = 0.98, 0.985, 0.99 and 0.995. 
The only significant deviation is at the lowest temperature 
T = O.OlJ/fcs. The fit function is given in Eqs. (||) with the 
parameters in Table ^ 

the data. The fit deviations from the 802 QMC and 
1749 TMRG data are shown separately in Figs. |2^(a) 
and |2^(b), respectively. A comparison of the two figm'es 
shows that the TMRG data are, on average, significantly 
more precise at a given temperature. 

After the parameters in the present x* (c^i t) fit function 
were finalized, as a check on the accuracy of the fit func- 
tion for a values close to the uniform limit, we carried out 
QMC x*{t) simulations for alternating-exchange chains 
of length L = 400 and 800, factors of four and eight longer 
than the chains for which QMC data were combined with 
TMRG data to determine the fit function, respectively. 
The simulations were carried out for a = 0.98, 0.985, 0.99, 
and 0.995 at temperatures 0.01 < t < A. Overall, the fit 
function was found to be in extremely good agreement 
with the QMC data. For OA < t < 4, the x*{a,t) fit 
function agreed with the simulation data to within about 
±5 X 10~^ or better. The deviations of the fit function 
from the data for 0.01 < t < 0.4 are shown in Fig. |2^, 
along with the error bars on the QMC data. As can be 
seen from the figure, the only significant deviation of the 
fit function from the QMC data in this t range is at the 
lowest temperature t — 0.01 for each of the four a 
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FIG. 23. (a) Magnetic susceptibility x versus tempera- 
ture T for the spin S — 1/2 antiferromagnetic alternat- 
ing-exchange Heisenberg chain with values of alternation pa- 
rameter 5 from to 1 as shown, where x is scaled by 1/J and 
T by J in contrast to Fig. The parameter J = Ji/(l-|-5) is 
the average of the two exchange constants Ji and J2 alternat- 
ing along the chain, (b) Expanded plot at low temperatures. 
These x* {S, t) plots were generated using our two-dimensional 
X* (a, t) fit function which was converted to the variables (5, t) 
using Eq. (16b). 



values. Because the fit deviations at this temperature 
remain upon increasing the length of the simulated chain 
from L ~ 400 to L — 800, these fit deviations are most 
likely due to inaccuracies in the fit function, as expected 
at this lowest fitted temperature. 

For compounds showing spin-Peierls or other types of 
second-order spin dimerization transitions, it is more ap- 
propriate to scale X by 1/J and T by J, where J is the 
average of Ji and J2, in which case the appropriate al- 
ternation parameter is 5 rather than a. It is straight- 
forward to convert our x*(a,i) fit funct ion to the form 
X*{5,t), where t = UbT/J, using Eq. (16b). We have 
done this and plot the x*('^>'^) fit function versus tem- 
perature for a series of 5 values in Fig. p3|(a). An ap- 
pealing monotonic progression of x* {5^ t) with increasing 
5 is seen in Fig. E3(a); an expanded plot at lower tem- 
peratures is shown in Fig. |2^(b) . This formulation of the 
fit function allows accurate estimates to be made of the 
temperature-dependent spin gap in compounds exhibit- 
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ing spin-dimerization transitions, provided that the near- 
est neighbor S = 1/2 AF alternating-exchange Heisen- 
berg model is appropriate to them. An illustration of 
the procedure and the results to be gained will be given 
later when we model the x(^) data for NaV205. 



VI. SPIN GAP FROM TMRG 



S = 1/2 AF Alternating-Exchange Chain 



According to Eq. (11c), if highly precise x*(0 data in 
the low-t limit are available, the spin gap A* can in prin- 
ciple be computed directly from the derivative of these 
data with respect to inverse temperature. However, in 
general the maximum temperature of the \ow-t limit re- 
gion is ill defined since it depends on how precise and 
accurate the data are and the accuracy to which A* is 
to be determined. Therefore, in practice one could de- 
fine a tem perature-dependent effective spin gap A*jj from 
Eq. (In]) as 



d{l/t) 



(59) 



and then try to ascertain A* from the extrapolated zero- 
temperature hmit A* = limt^o A*fj(i). Using Eq. (lit) 
would be less desirable and precise because a fit of this 
type typically averages A*jj(i) over a rather large tem- 
perature range. 

An overview of A*g(a, t) determined from our TMRG 
X*(a, t) data for 0.8 < a < 0.995 using Eq. (!§) is shown 
in Fig. p^(a). At the lowest temperatures, and for a 
not too close to 1, the A*g(a, t) data do approach a con- 
stant value with decreasing confirming the applicability 
of Eqs. ([l0| ) and prior assumptions and hence Eqs. ( |ll| ) 
and ( p9[ ) to the alternating chain, and the approximate 
values of A* (a) can be estimated from the figure. Closer 
inspection reveals that A*g(a, t) shows a weak maximum 
before decreasing by ~ ^% to A* as t ^ 0, as illustrated 
in Fig. ^(b) for a — 0.8. For this among other rea- 
sons, we will not use Eq. ( p9| ) to extract the spin gaps 
from our TMRG x*(c*iO data. On the other hand, we 
need to know whether such behavior is expected, since 
it could conceivably arise from systematic errors in the 
TMRG calculations. Therefore, in the next section we 
study the A*g(a, t) expected at low temperatures for the 
alternating-exchange chain. As part of this study, we for- 
mulate and discuss the fit function which we will use in 



Sec. |VIB| to extract A* (a) from our TMRG %*(«, t) data 
at low temperatures. 



A. Effective spin gap A*jj(A*,t) for the 
alternating-exchange chain 

From our definition of A*ff in Eq. ( ^9[) , a discussion of 
how this quantity varies with t at low t requires an inde- 
pendent estimate of x*(q^, i) for the alternating-exchange 
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FIG. 24. (a) Overview of the effective spin gap /S.*g(a,t) 
vs temperature t for the S = 1/2 AF alternating-exchange 
Heisenberg chain, derived from our TMRG x*(q^i *) data using 
the definition in Eq. (jH^), where a = J2/J1 is the alternation 
parameter, (b) Expanded plot of A*g (t) for a = 0.8 at low 
temperatures from (a). 

chain, which must include at least the leading order cor- 
rection to the low-t limit in Eqs. (|l0|). As a first attempt, 
we used the general expression for x*(0 Eqs. (^, which 
requires as input the one-magnon dispersion relation e(fc) 
for the alternating chain. For this we used the explicit 
£(A*,/c) for < a < 1 in Eqs. (||) that we presented 
and discussed previously. The resultant x*(^) is plotted 
for eleven A* values in Fig. |2^(a), where the results are 
designated by X*(^) in the figure. Although x*^'^K^*,t) 
is exact in both the low- and high-t limits, the results 
are only qualitatively correct at intermediate tempera- 
tures, as can be seen by comparing Fig . |2^(a) with the 
QMC and TMRp.data and fit in Fig. |o[ Troyer, Tsunet- 
sugu and WiirtzEII obtained a very good fit oix*^^K^*,t) 
to QMC X* {i) simulation data over a large temperature 
range for the S = 1/2 two- leg Heisenberg ladder with 
spatially isotropic exchange; however, they assumed a 
e(A*, k) in the fit function which was later found to be 
inaccurate over much of the Brillouin zone. 

We formulated an approximation [designated as x*^^''] 
which is more accurate in the low-temperature range, and 
which we will use in the next section as a fit function to 
fit our TMRG x*(i) data at low t to extract A* (a). The 
function x*^^"* was obtained by summing the susceptibil- 
ities of isolated dimers with a distribution of singlet-trip- 
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FIG. 25. Magnetic susceptibilities x'^' (a) and x'^' (b) vs 
temperature T for the spin 5 = 1/2 antiferromagnetic alter- 
nating-exchange Heisenberg chain, calculated using two dif- 
ferent approximations for X*(^)i respectively (see text). Note 
the different scales for the ordinates in (a) and (b). 

let energy gaps given by our one-parameter dispersion 
relation e(A*,A:) for < A* < 1 in Eqs. (|2|), which 
takes into account the interdimer interactions. Thus from 
Eq. (pa|) we simply obtain 



*(2) 



dk 



3 -|_ ee(A*,fc)/t 



(60) 



Note that we make no assumptions here about the form 
of the function A* (a), since only A* appears in the ex- 
pression. This x*^^)(A*,i) is exact in both the low- 
and high-i limits, as is both repro- 

duce x*{i) for the isolated dimer (A* — 1) exactly, but 
;!^*(^)(A*, t) is more accurate at intermediate tempera- 
tures for a < 1 as shown in Fig. |2^(b). In addition, by 
comparing and x*^^^^*.^) with the TMRG 

X*(a,t) calculations at low we found that the low-i 
corrections to the low-i limit in Eqs. ( p^ ) are much more 
accurately given by X*^'^Kt) than by x*''^\t). We wiU 
therefore not discuss x*''^\i) further here. 

At low temperatures, the approximation x*''^'' (A*, i) is 
expected to accurately describe the leading-order t cor- 
rections to the \ovf-t limit only as long as the average 
number of magnons rim occupying a state near the min- 
imum in the one-magnon band is much less than unity. 
Using the expression for the boson occupation number 
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FIG. 26. Effective spin gap Acff, defined in Eq. (|5S|) 
and computed using Eq. (pOl), vs reduced temperature 
k^T / A = f/A* for the S = 1/2 antiferromagnetic al- 
ternating-exchange Heisenberg chain for various values of 
A* = A/Ji, where t = k^T/Ji. A limiting behavior is seen 
for A* 0, for which the maximum in Aeg (fcer/A) occurs at 
A:bT/A ^ 0.14. The hnear in T behavior of Acft (fcBT/A 0) 
for the isolated dimer (A* = 1) is due to the identically zero 
width of the one-magnon dispersion relation for this A* value. 



for this case. 



gA-/t _ 1 



(61) 



yields t/A* = 0.22 and 0.42 for 0.01 and 0.1, re- 

spectively. Thus, when fitting our \ow-t TMRG 
data by the fit function x*^^H^*; the following sec- 
tion, we expect x*(^)(A*,i) to be sufficiently accurate 
only for t/A* < 0.4. For this reason, our fits will be 
limited to this maximum scaled temperature. 

We have computed Aeff(A*,T)/A from x*'^^^*.*) in 
Eq. (|60|), using the definition in Eq. (p9|), and plot the 
results vs fceT/A in Fig. |2|. For the dimer (A* = 1), 
one finds analytically that Aoff(T)/A = 1 - 2kBT/A to 
lowest order in T. On the other hand, for < A* < 1, 
the initial dependence is positive and quadratic in T, and 
a maximum is seen in Aoff(T)/A, which for A* < 0.4 
occurs at t/A* = ksT/A « 0.14 with a height of w 0.5%. 
This height is quantitatively consistent with the data in 
Fig. |2|(b) derived from the TMRG for a = 0.8. 

Thus the weak maximum seen in that figure is not a 
spurious effect. 



B. Fits to the low-t TMRG x*{ct,t) data 

We were tempted to fit A*^^ (a, t) derived from the low- 
t TMRG x*(Q!j^) data, as discussed above, to obtain 
the spin gaps A*{a,t). However, this procedure would 
have weighted the x*(«ji) data in an ill advised way. 
We therefore decided to do conventional fits of the low-t 
X*{a,t) data by the fit function x*^^H^*,t) in Eq. (|6^). 
For a given a, this is a one-parameter (A*) fit function 
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and the fits are therefore stringent tests of both the ap- 
propriateness of the fit function and the precision and 
accuracy of the data. Because the temperature depen- 
dence of the accuracy of the calculations is unknown ex- 
cept for the uniform chain data (see Fig. |l^), we assumed 
that all data for a given a in a given fitted temperature 
range have the same accuracy. Thus in the nonlinear 
least-squares fits for each a we minimized the square of 
the relative rms deviation of the fit from the data 
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1 



E 



(62) 



where Np is the number of data points fitted, which was 
usually between 250 and 1500. 

Due to the presence of the spin gap A* in the expo- 
nential of the fit function, tJrms is extremely sensitive to 
the precise value of A* when low-t fits are carried out. 
For example, close to the optimum A* fit value, a change 
in A* by only 0.0001 0.1 %) can change arms by up to 
~ 300%. Thus a few percent accuracy in the x*(i) data 
at low t is sufficient to allow A* to be determined for a 
given fit to a precision better than 0.0001. For a given 
a, the obtained A* was found to be insensitive, typically 
to within « 0.0002, to the t range of the fit, as long as 
the maximum fitted temperature satisfied t/A* < 0.4, 
consistent with the above discussion of the boson occu- 
pation number. This lack of sensitivity of the value of the 
fitted A* to the precise fit range demonstrated that the 
fit function x*^^^ (t) is an appropriate one. Depending on 
the a value and the t range of the fit, typically 
between 0.1 % and several percent. 

The A* (a) values obtained from the fits are listed 
in Table 
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together with the estimated accuracies in 
parentheses. Note that a quoted accuracy is associated 
with variations in A* in fits to a specific set of data for 
a given a over various temperature ranges as discussed 
above, and does not include possible systematic errors 
due to, e.g., the finite fixed number of states kep t in the 
TMRG calci|aa|tiMif>_ Also included in Table |l| are lit- 
erature dataHtZIH'Ld which will be compared with the 
present results in the next section. Log-log plots of the 
low-t data and fits are shown in Fig. where on the 
scale of this figure, for most a values the data and fit 
are identical (cannot be distinguished) within the fitted 
temperature range. Extrapolations of the fits to higher 
and lower temperatures are also shown for comparison 
with the data. 



VII. COMPARISONS OF THE CALCULATIONS 
WITH PREVIOUS WORK 

A. Spin gap 

Our A* ((5) spin gap data determined by fitting our 
TMRG X (t) data by Eq. @) are plotted in Fig. ||(a). 




FIG. 27. Log-log plots of reduced magnetic susceptibil- 
ity X* = xJ^/^g^f^B "^s reduced temperature t = ksT/Ji 
(solid circles) for spin S = 1/2 antiferromagnetic alter- 
nating-exchange Heisenberg chains with alternation param- 
eters a = J2/J1 shown in the figure, calculated using 
TMRG. The corresponding fits to the lower temperature data 
[kBT/A{a) < 0.4] by Eq. I p(\ ) are shown as the dashed curves, 
which are extrapolated to lower and higher temperatures in 
the figure. The discontinuity in the data and fit for a — 0.995 
at t « 0.0043 is due to an increase with increasing t at that t 
in the number of states kept in the calculations from 150 to 
256. The spin gaps A* (a) = A{a)/Ji found from the fits are 
given in Table pil 

along with the results of previous workersElEBil'S listed 
in Table |ll^. The solid curve is the function A* = 25^^^ 
in Eq. r ^Obj) proposed by Barnes, Riera, and Tennant 
(BRT).EjThe overall behavior of the data in Fig. |2|(a) is 
well described by this function, but significant deviations 
of the data from the curve occur as illustrated in the 
expanded plot for S < 0.1 in Fig. ^(b). The error bars 
are included with each plotted data point in Fig. |28|(b), 
except for the data of Ref. ^ which were not available, 
but they are all small and not clearly seen. Our values 
forjL< 0.1 are significantly smaller than those of Uhrig et 
a/.J£3 where the differences are far outside the combined 
limits of error, and are larger than tho se of Au gier et a/.ES 
As will be seen explicitly in Sec. VIII C below, our 
X*{5,t) fit function allows S{T) to be determined for real 
materials by using the fit function to model experimental 
x{T) data. However, if one would like to determine the 
spin gap A(T) from the derived S{T), an expression is 
needed for A(d) over the entire range < (5 < 1 in order 
to be generally useful and applicable. At present, the 
only extant expression is that of BRT in Eqs. (|l9|). As 
seen in Fig. ^ and in Table IV below, this expression 
is only an approximation that fits neither BRTs' A* (a) 
data for 0.1 < a < 0.9 nor our TMRG spin gap data for 
0.8 < a < 0.995 to within the respective error bars. To 
formulate a more flexible expression, we modify BRTs' 
formula to read 



A*{S) = 



m 
J 



(63a) 



so the ^-dependent power y is 
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TABLE III. Spin gaps A* (5) = A(5)/J and A* (a) = A(a)/Ji for the 5* = 1/2 antiferromagnetic alternating-exchange 
Heisenberg chain as determined using T = DMRG calculations by Uhrig et al. [A*u{5) and A*u(a)1 (Ref. |6^ , by Barnes, 
Riera, and Tennant [A*brt(5)] using multiprecision methods (Ref. by Augier et al. [A*a((5)] (Ref. and by us [A*pw(5) 
and A*pw(a)] in the present work (pw) from our TMRG x*(«ii) data as described in the text. Two A values are given 
for a — 0.995 in the present work: the first (larger) value is for the number of states kept in the calculations m — 150 at 
t < 0.004533, whereas the second (smaller) value is for m = 256 at t > 0.004533. The A* (a) data of Barnes, Riera, and 
Tennant are given in Table IV, Additional literature data include those of Ladavac et al. obtained using a Green's function 
Monte Carlo technique on rings of 6 to 200 spins (Ref. I^): A^{S) = 0.1815(5) {S = 0.04), 0.2156(1) (0.05), 0.301(1) (0.08) and 
0.3603(1) (0.10). 
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FIG. 28. (a) Spin gap A/ J vs alternation parameter 5 for 
the S = 1/2 antiferromagnetic alternating-exchange Heisen- 
berg chain. Our data (•) were determined by fitting our 
TMRG X* (*) data by Eq. mu) and are shown along with data 
of Barnes, Riera, and Tennant (Ref. (filled squares), Uhrig 
et al. (Ref. (o), Ladavac et al. (Ref. (open squares), 
and Augier et al. (Ref. |58[) (open diamonds). The solid curve 
is the function (Ref. |jA/J = 2S^^*. (b) Expanded plot of 
the data and curve in (a) for 5 < 0.1. Error bars for the data 
are not shown in (a), but are shown in (b) for all data except 
for those of Augier et al. 



ln[A((5)/2J] 
"2" in Eq. 



(63b) 



The numerical prefactor "2" in Eq. (63a) must be re- 
tained in order to reproduce the exact A*{S — 1) — 2. 
Shown in Fig. p9|(a) is a semilog plot of y versus 6 for 
the same numerical data as in Fig. |2^. This plot [and 
Fig. ^9|(b) below] explicitly shows, from BRTs' data, 
that the exponent deviates significantly from the value 
3/4 even for S ^ 1. The plot also clearly differentiates 
the various numerical data for small 6 by the different 
groups, and shows that one of our two data points from 
the TMRG for S = 0.0025 (the one derived from m = 256 
data at high t) is not in agreement with the trend of the 
remainder of our data. This data point will not be in- 
cluded in the plot and fit to be discussed in the next 
paragraph. 

Our y{S) data at small 6 are in agreement with both 
the magnitude and trend of BRTs' data at larger S. The 
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FIG. 29. (a) Semilog plot vs alternation parameter S of the 
exponent y = ln(A/2 J)/ In (5 in the expression A/ J = 2S^ 
for the spin gap A of the 5 = 1/2 antiferromagnetic alter- 
nating-exchange Heisenberg chain. The data and symbol ref- 
erences are the same as in Fig. Each data point has an 
attached error bar except for those of Augier et al. (Ref. [s^ ) 
(open diamonds), (b) Expanded view of our y{S) data (•) and 
those derived from the numerical spin gap data of Barnes, Ri- 
era, and Tennant (Ref. (filled squares), along with the fit 
in Eqs. (33c) and (33d) to the two combined y{5) data sets 
(solid curve). 



y{S) for these two data sets from Fig. ^(a), with the 
exception of one of our two data points for S — 0.0025 just 
noted above, are plotted together on an expanded vertical 
scale in Fig. ^(b) where a rather smooth behavior of y{6) 
is seen over the combined range of the two calculations 
0.0025 <S<1. With the behavior in Fig. |2|(b) in mind, 
we formulated a five-parameter fit function for these two 
combined y{S) data sets that yields the correct limits 
A* ((5 ^ 0) = and A*{S ^ 1) = 2, with the property 
lima^o y{S) — const, given by 



y((5) — y{l) + ni tanh 



n2 tanh 



mi V 1712 ) 

i^lnfi^) 
TOi V mi I 



(63c) 



An unweighted fit of this expression to all the data in 
Fig. |2^(b) yielded the parameters 
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FIG. 30. Spin gap A/J vs 5^/^! In^l^^^ obtained from the 
data in Fig. Eq. The straight hne passes through the origin 
with slope 3.3. 

2/(1) = 0.74922 , m = 0.00776 , 7^2 = -0.00685 , 



mi = 3.3297 



m2 



-2.2114 



so that 



lim y{S) = y(l) - ?ii + = 0.7346 

(5^0 



(63d) 



(64) 



The fit is plotted as the sohd curve in Fig. ^(b). As 
can be seen from the figure, our data are fitted to within 
our error bars. I n ad dition, when the y{6) fit f unction in 
Eqs. (63c) and ( |63d|) is inserted into Eq. (63a), the pre- 
dicted values of A* ((5) are in agreement with each of the 
values of BRT at larger S to within 0.0001, which is suffi- 
cient for modeling experimental data. The 6 = limit of 
y{S) in Eq. ( |64| ) is in agreement with the theoretical ef- 
fective value y(0) = 0.72(3), which was obtained atithout 
the log correction term by Singh and WeihongO from 
an eleventh-order dimer series expansion of the triplet 
dispersion relation. We will use Eqs. ( |6^ ) to compute 
A(r) f rom the experimentally derived SJT) for NaV205 
in Sec. 
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In order to test the critical behavior prediction A* = 
A(52/3/| ln<5|i/2 in Eq. (|2^), which need only hold in the 
asymptotic critical regime (5 in contrast to the fit 
function for < (5 < 1 in Eqs. (|6|), in Fig. |^ is plot- 
ted A/J vs .52/3/1 ln5|i/2 in the region 6 < 0.06 for the 
same data and symbols as in Fig. A proportional- 
ity appears to be developing in our data for S ^ 0.005, 
as shown by the straight line with slope A = 3.3 pass- 
ing through the origin of the figure, suggesting that the 
asymptotic critical regime begins with decreasing 6 be- 
low S - 0.005 (a > 0.99). High-accuracy A*{S) data 
for S < 0.001 are needed to test this conjecture. From 
Fig. the slope 3.3 of the line drawn is evidently a 
lower limit of the prefactor A within the actual asymp- 
totic critical regime. 
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FIG. 31. (a) Magnetic susceptibility x versus tempera- 
ture T for the spin S — 1/2 antiferromagnetic alternat- 
ing-exchange Heisenberg chain with alternation parameter 
a = 0.2, 0.4, 0.6, 0.7 and 0.8 (symbols), calculated by Barnes 
and Riera (Ref. y) using the Lanczos technique. Our x*(i) fit 
function as in Fig. 1^ (soiid curves) for the same a values is 
shown for comparison, (b) Deviation of the data of Barnes 
and Riera from our fit function vs T. 



B. Numerical x*(Q:,t) results 

Barnes and Riera previously carried out exact diag- 
onalizations of Hamiltonian (|lj) for 5=1/2 alternat- 
ing chains-|Of length up to 16 spins using the Lanczos 
techniqueJa Their computed x*(/) values for a = 0.2, 
0.4, 0.6, 0.7, and 0.8 were extrapolated to the bulk limit 
and the results are shown as the symbols in Fig. pT|(a). 
Our fit function as in Fig. 20 for the same a values 
is plotted as the solid curves in Fig. |3^(a), which are 
seen to be in good overall agreement with the calcula- 
tions of Barnes and Riera. The deviations of the data 
of Barnes and Riera from our fit function are plotted 
vs temperature in Fig. ^(b). The average deviation of 
their data from our fit function is very small for each data 
set: -0.41, -^0.33, -0.40, -0.26, and -1-0.79 x 10"^ for 
a = 0.2, 0.4, 0.6, 0.7, and 0.8, respectively. The absolute 
rms deviations Crms of their data from our fit function for 
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a = 0.2, 0.4, 0.6, 0.7, and 0.8 are (in units of IQ-^) 1.73, 
1.43, 0.73, 0.78, and 3.76, respectively. We conclude that 
their data are in good quantitative agreement with our 
data and fit, with the exception of their data point for 
a = 0.8 at their lowest temperature t = 0.05. 



C. Bulaevskii Theory 

BulaevskiiB calculated analytically in the 

Hartree-Fock approximation. He first obtained an in- 
tegral equation for the magnon spectrum E{k): 



e{k) 
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Ji 



■\/l + — 2a cos k 
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(65) 
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(66) 



where k is measured in units of 2'k j d. d = 1 is the lattice 
repeat distance along the chain, which is twice the aver- 
age distance between spins. He then expressed x*(^) i'^ 
terms of e(fc): 



m ^ ^ r 



2 + {l + a)F{t) ' 
dk 



(67) 



/o cosh^[e(fc)/(2t)] ■ 
At i = and a ^ 0, from Eqs. (|6q) we obtain 



Ci(a) = i<j (1 + q)E 



4q 



(1 - a)K 



(l + a)2 
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1 



C2(a) = -<^ (l + a)E 



(1 - a)K 
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(68) 



(l + a)2 

4Q! 



(l + a)2 



where K(j/) and E(?/) are, respectively, the complete el- 
liptic integrals of the first and second kinds. The disper- 
sion relations versus a at t = are obtained by inserting 
Eqs. (pq ) into (65) and a selection of results is shown in 
Fig. ppFrom Eqs. (||) and (||), at t = the spin-gap 
Afc^oTo^) at fc = is given by 



LU 




FIG. 32. Dispersion relations E{k) at temperature T = 
from Bulaevskii's theory (Ref. |^ in Eqs. ( |65| ) and (^) for ten 
values of the alternation parameter a. The inset shows an 
expanded plot near = of E{k) for a = 0.9. 
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(69) 



This expression gives the actual spin-gap for < a < 
0.79. However, for 0.79 ^ a < 1, the minimum in the 
dispersion relation does not occur at A: = 0, as illustrated 
in an expanded plot of E(k) for a = 0.9 in the inset to 
Fig. The wave vector ko at which the minimum spin 
gap Ab occurs is plotted versus a in Fig. |3^(b). The Ab 
from Bulaevskii's theory at i = is plotted versus a as 
the solid curve in Fig. 33(a), and a few representative val- 
ues are given in Table I^. The predictions of Bulaevskii's 
theory are in very eopd agreement with those of Barnes, 
Riera, and Tennanta for a < 0.3, but the agreement be- 
comes progressively worse as a increases further. From 
Eqs. ( pSi ) and (|6^), E{k) is temperature dependent. In 
addition, in the range 0.79 < a < 1 for which kQ ^ Q 
at i = 0, we find that fee depends on t, as shown in 
Fig. ||. From Fig. ||, /cq ^ at t « 0.083, 0.122, 0.131, 
0.125, and 0.111 for a ^ 0.8, 0.85, 0.9, 0.95, and 0.99, 
respectively. 

We computed x*it) by inserting e{k) in Eq. ( |65| ) into 
Eqs. (|66|), numerically solving the latter two simultaneous 
equations for C'l and C2 at each t, and then inserting the 
resulting e{k) into Eqs. (|67|). The progression of x*(t) 
with increasing a from ft-OOl to 0.99 is shown in Fig. pB. 
As noted by BulaevskiiJH the values of x* at the maxima 
are too large and the temperatures at which these occur 
are too small by ~ 5-10% (compare Fig. ^with Fig. po| ). 

At low temperatures 0.033 < i < 1/4, Bulaevskii fitted 
X*(a, t), calculated from Eqs. (|67|), by the two-parameter 
form 
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FIG. 33. (a) Energy gap A versus alternation parameter a 
for the 5 = 1/2 alternating chain, as calculated by Barnes, 
Riera, and Tennant (Ref. p3|) (•) and by us using the theory of 
Bulaevskii (Ref. ^) (solid cu rve) . The dashed curve is a plot 
of A versus a given in Eq. (19a). The values of A obtained 
by Bulaevskii (Ref. ^) by fitting his numerical calculations of 
X(r) for 0.033 < fcBT/Ji < 0.25 according to Eq. (0) are 
shown as the open squares, (b) Wavevector fee, at which the 
minimum spin gap occurs in the magnon dispersion relation 
at T = 0, vs alternation parameter a. 
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FIG. 34. The temperature T dependence of the wavevec- 
tor ko at which the spin-gap occurs in the triplet magnon 
dispersion relation (^), for five values of the alternation pa- 
rameter a. 
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FIG. 35. Magnetic susceptibility x versus temperature T 
for 0.001 < a < 0.99 as predicted by the theory of Bulaevskii 
(Ref. |) in Eqs. (||)-(|3). 
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(70) 



and obtained values of ae and AB,Fit/^i for < a < 0.9 
which are reproduced in Table [V; AB,Fit(Q!)/Ji is plot- 
ted as the open squares in Fig. |33|(a). Note that the tem- 
perature exponent in the prefactor to the e xponential is 
7 = 1, contrary to the 7 = 1/2 in Eq. ( |lOb| ) which is ex- 
pected in the low-t limit for any ID S' = 1/2 Heisenberg 
spin system with a spin gap (and with a nondegener- 
ate one- magnon band with a parabolic minimum). We 
have confirmed that over the temperature range fitted by 
Bulaevskii, one indeed obtains 7 ~ 1 for the best fit of 
Eq. (70) to numerical calculations of x*(a,t). We infer 
that the discrepancy between Bulaevskii's 7=1 and the 
expected 7 = 1/2 arises because the fits were not car- 
ried out completely within the low-i limit. This issue is 
discussed in more detail below. 

Equation (|7^) together with Bulaevskii's table of 
{03(0;), AB,Fit(Q;)/Ji} values were subsequently used ex- 
tensively in the analysis of experimental x{T) data for 
compounds exhibiting spin-Peierls transitions to deter- 
mine the alternation parameter a at low temperatures 
T <ti Tc where the experimental spin-gap is nearly inde- 
pendent of T. However, from Tajble IV and Fig. ^(a), 
the AB,Fit(Q:) values of BulaevskiiQ are in generally poor 
agreement with the actual spin gaps Ab(q;) of his the- 
ory and with those [AbrtIck)] calculated for the same a 
values by Barnes, Riera, and TennantB Therefore, one 
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should consider the Ab^fh parameters as fitting param- 
eters only, with no direct relation to the actual spin gap. 

According to Eq. (|8^) for x* (^) of the isolated dimer 
which is a zero-dimensional spin system, the form of 
X*{i) in Eq. ( |70| ) with 7 = 1 is correct for a = and 
t ^ 0. On the other hand, for one-dimensional spin sys- 
tems such as the two-leg spin ladder (and the alternating- 
exchange chain) at temperatures /cbT <C A and k^T <ti 
one-magnon bandwidth, Eqs. ( p^ ) apply, with 7 = 1/2, 
assuming that the triplet one-magnon dispersion rela- 
tion E(k) is parabolic at the minimum. In this case one 
expects 7 = 1/2 at sufficiently low t for any finite a. 
Thus, in the temperature region of v alidit y of Eq. ( |10aD , 
a plot of the left-hand-side of Eq. (11a) vs Int should 
give a straight line with slope —7. Shown in Fig. [s^ are 
such plots, obtained using our x*(Q;,i) calculated from 
Bulaevskii's theory as described above, for a = 0.001 
to 0.99. For a = 0.001, a crossover is clearly evident from 
7 = 1 to 7 = 1/2 with decreasing t. The other curves also 
exhibit signs of a crossover, with 7 w 1/2 at the lowest 
temperatures, with the exception of the curve for a = 0.8. 
For this a value, which is just above the value a ~ 0.79 at 
which kQ becomes nonzero at t = [see Fig. |3^(b)] , the 7 
at the lowest t is intermediate between the values of 1/2 
and 1, and the assumption of a parabolic form for E(k) at 
the band minimum is evidently not satisfied (see Fig. 32). 
In fact, Troyer, Tsunetsugu and Wiirtztil calculated the 
low-i limit of x*(i) for ID systems with general dispersion 
relation e(k) = A* -|- c*|A:a|"', where k is the deviation of 
the wave vector from that at the band minimum. They 
found the same form x*(t) = (An/f) exp(— A*/i) as for 
the parabolic case n = 2, but where 7 = 1 — (1/n). Thus, 
e.g., 7 = 2/3 and 3/4 for n = 3 and 4, respectively. This 
range of 7 is consistent with the slope of the data at the 
lowest temperatures for a = 0.8 in Fig. p6[ 

The predictions of Bulaevskii's theory for x*(i) from 
Fig. ^ are compared with our x*(a, t) fit function (solid 




FIG. 36. Log-log plot of x*(i) e^'^''''' at low t ver- 
sus redu ced temperature t for 0.001 < a < 0.99 [see 
Eq. ( |lla|) 1 as p redicted by the theory of Bulaevskii (Ref. ^) in 
Eqs.l|)-(|5. 



curves) for a = 0.2, 0.4, 0.6, 0.8 and 0.99 (as in Fig. gO|) in 
Fig. p^, where the Bulaevskii prediction for each of these 
a values is shown as the corresponding dashed curve. 
The disagreement between the two calculations becomes 
progressively more severe as temperature decreases and 
as the uniform chain limit is approached with increas- 
ing a. Therefore, the accuracies of the a and Ji values 
previously extracted from experimental data at low T 
for compounds with a ^ 1 using Bulaevskii's theory are 
unclear. Our x*(Q!,t) fit function now provides a much 
more accurate and reliable means of extracting exchange 
constants and spin gaps from experimental x(r) data. 



TABLE IV. Prefactor a and spin-gap A describing the low-temperature spin susceptibility in Eq. ( |7C| ) of the S = 1/2 
alternating chain with alternation parameter a [and S = {1 — q)/(1 +cy), see Eqs. (p^J. Fit parameters given by Bulaevskii 
(Ref. ^ (ob, Ae.Fit) that he obtained by fitting to his \ow-t x*(*) calculations using Eq 
actual spin gap values Ab/Ji in Bulaevskii's theory by numerically solving Eq s. (m) and (IM) 



rOD are shown. We obtained the 
Also included are the accurate 

calculations of the spin gap Abrt/Ji by Barnes, Riera, and Tenn ant (Ref. |53[), which are compared with numerical values of 



their approximate form A(a)/Ji « (1 - ay'^{l + ay''^ [Eq. (|l9aD 
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FIG. 37. (a) Magnetic susceptibility x versus temperature 
T for alternating chains with a = 0.2, 0.4, 0.6, 0.8, 0.9, and 
0.99 (solid curves) generated using our x*(ai *) fit function as 
in Fig. These are compared with the predictions of the 
theory of Bulaevskii (Ref. ^ (corresponding dashed curves), 
(b) Expanded plots at low T from (a). 



VIII. MAGNETIC SUSCEPTIBILITY OF NaVaOs 

Crystals of Nag ggg(3)V205.oo(6) were grown at the 
Max-Planck-Institut fiir Festkorperforschung, Stuttgart, 
in a Pt crucible in flowing Ar atmosphere by a self- 
fliox rnethod from a 5:1:1 mixture of NaVOs, V2O3 and 
V205.La The flux was dissolved by boiling the solidified 
melt in distilled water. X-ray powder diffraction pat- 
terns collected with a STOE diffractometer yielded the 
lattice parameters a = 11.3187(8) A, b = 3.6111(3) A, 
and c = 4.8007(5) A. Chemical analyses on two inde- 
pendent representative samples of the batch were per- 
formed with a standard AAS analysis technique for V 
and Pt and ICP emission spectroscopy for the Na con- 
tent. The oxygen content was determined by measuring 
with IR spectroscopy the amount of CO generated when 
the sample is fused in a graphite crucible at 2700 °C in 
vacuo. Platinum impurities above the level of sensitivity 
of the analysis (500 ppm with respect to V) could not be 
detected. 




NaV^Og Crystals 
H = 1 T parallel to V chains 



) 100 200 300 400 500 600 700 800 
T(K) 

FIG. 38. Magnetic susceptibility x in 8- field H = IT par- 
allel to the V chains {b axis) versus temperature T for two 
crystals of NaVaOs as indicated. The solid curve is a "Fit" of 
the data by the theoretical prediction for the S — 1/2 uniform 
Heisenberg chain with parameters in Eq. (fz^). 



At Ames Laboratory, single pjystals of NaV205 were 
grown out of the ternary meltJZJ Powders of V2O5 and 
V2O3 were prepared by oxidizing and reducing NH4VO3 
at 600 °G and 900 °C, respectively. The resulting V2O5 is 
reacted with Na2C03 at 550 °C yeflding NaV03. About 
10 grams of NaV03, V2O5, and V2O3 in the molar ra- 
tio 32:1:1 were placed in a Pt crucible and sealed in an 
evacuated quartz tube. The melt was then slowly cooled 
from 800 to 660 °C over 50 h and the remaining liquid 
was decanted. Small amounts of solidified melt remaining 
on the crystals were dissolved with hot water. Typical 
dimensions of the ribbon-shaped crystals grown in this 
manner are 0.5 x 1.5 x 11 mm'^ with the c axis perpendic- 
ular to the plane of the ribbon, the b axis along the length 
of the ribbon and the a axis along the width of the ribbon, 
with lattice parameters a ~ 11.303 A, b « 3.611 A, and 
c « 4.752 A. The crystal denoted as ALl has a mass of 
8 . 2 mg and approximate planar dimensions 1.5x2.5 mm^ . 

The magnetic susceptibflity xiT) = M{T)/H of the 
crystals was measured using Quantum Design SQUID 
magnetometers at Stuttgart and Ames. The measure- 
ments on eight crystals of NaV205 in Stuttgart were car- 
ried out in a field H = IT along the V ladder (6) axis 
direction in various temperature ranges between 2 and 
750 K. Measurements of the anisotropy of x(^) along 
the a, 6, and c axis directions were carried out from 2 to 
300 K in Ames on crystal ALl in = 2T. 

The results for two of the crystals up to 750 K are 
shown in Fig. The data illustrate the variabilities 
we have observed between measurements along the same 
axis on different crystals. Above ^ 50 K, the two data 
sets are nearly parallel, with the difference between them 
being w 3 — 4 x 10^^ cm"^/mol; we have no explanation for 
this difference, and no comments have been made in the 
literature about such variabilities and/or their origins in 
x{T) along the same axis in different crystals that we are 
aware of. The data from Tc w 33-34 K up to 300 K are 
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in approximate agreement with the single crystal data of 
Isobe, Kagumi and Ueda taken in this T range along the 
same axis in _ff = 5 T.O A variable Curie- Weiss-like con- 
tribution x^^(T) to x(r) occurs below - 20K which is 
attributed to paramagnetic defects, impurities, inclusions 
and/or intergrowths in the crystals. Th e "Fit" shown in 



the figure will be discussed later in Sec. VIII B, 

The experimental data are analyzed with the general 
expression 



x(T) = xo + x^^(r) + x^p'"(T) 



^corc_^^VV 



imp 



T-9 ' 



(71a) 



(71b) 



(71c) 



(71d) 



where xo is the sum of a temperature independent 
and (nearly) isotropic orbital diamagnetic core contri- 
bution and a usually anisotropic and temperature inde- 
pendent orbital paramagnetic Van Vleck contribution. 
We estimate x™™ using the values —5, —7, —4, and 
-12 X 10-*5r«nVmol for Na+i, V+^ and 0"^ 

respectively,E3 yielding the isotropic value 



X 



-7.8 X 10" 



mol NaV205 



(71e) 



The second term in Eq. (71a) is the above- noted Curie- 
Weiss impurity and/or defect contribution and the last 
term is the intrinsic spin susceptibility, each of which may 
or may not be anisotropic. For a Heisenberg spin sys- 
tem, X* is isotropic, and therefore so is yf^™ apart from 
anisotropy in the g factor. The impurity Curie- Weiss 
term x^^ [T) can be anisotropic if the impurities are de- 
fects or intergrowths in the crystals which have atomic 
coordination principal axes which are fixed with respect 
to the crystal axes rather than being randomly oriented. 
We model our x.{T) dat a according to Eq. ( [71a[ ) in terms 
of the x*(i) in Eq. (71d), which are (fit functions to) the- 
oretical susceptibility calculations presented in previous 
sections. Before moving on to do that, we first experi- 
mentally examine the anisotropy in x{T) of NaV205 and 
its implications in the next section. 



A. Anisotropy of the magnetic susceptibility 

The magnetic susceptibilities of NaV205 crystal ALl 
along the a, 6, and c axes are plotted vs temperature 
in Fig. |39|(a), where the a and c axes are perpendicu- 
lar to the V chains which run along the h axis, and the 
c axis is perpendicular to the trellis layers that the V 
chain/ladders reside in. The data are similar to the aniso- 
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FIG. 39. (a) Magnetic susceptibility x versus temperature 
T in a field H = 2T parallel {H\\b) and perpendicular {H\\a, 
H\\c) to the V chains in NaV205 crystal ALl. (b) The data 
in (a) corrected for the respective Curie- Weiss contributions 
^cw _ (j^ _ attributed to paramagnetic defects or 

impurities. 

tropic x{T) data reported by Isobe, Kagami and Ueda,El 
although the anisotropics we measure at both room tem- 
perature and at low temperatures are somewhat larger 
than they reported. The anisotropies at low tempera- 
tures are seen more clearly if the respective impurity term 
X^^(T) in Eq . (|71cD is subtracted from each data set, as 
shown in Fig. |39|(b). The impurity Curie constant Cimp 
and Weiss temperature 9 for each direction of the applied 
field were determined by the requirement that xiT) be- 
come independent of T for T ^ 0. The fitted values of 
Cimp were found to be slightly anisotropic and are given 



in Table Vl below. The values of Cimp are equivalent to 



the contribution of only 0.07 mol% of S* = 1/2 impurities 
with g — 2; if the impurity spin is actually greater than 
1/2, the concentration of paramagnetic impurities could 
be much less than this estimate. From a comparison of 
Figs. |9|(a) and ^(b), x^^iT) is seen to make a negli- 
gible contribution to the measured xiT) above ~ 100 K. 
Since in the presence of a spin gap x^^^'" = at the lowest 
temperatures, from Eqs. ([7l| ) and Fig. ^(b) we obtain 
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Xb 



Xa 



: 18.7 



X 10 ^ 

mol 



: 20.0 X 10'^ ^ 
mol 



xr - 13.3 
(T « Tc) . 



X 10"^ 

mol 



(72) 



From Fig. |3g|(b), the anisotropies of xiT) are seen to 
be quite temperature-dependent upon heating through 
Tc = 33.4 K. These results are surprising, because ;y^p'" is 
expected to be isotropic (apart from the small anisotropy 
due to the anisotropic g factor), with x''^'"(T' 0) = 
because of the spin gap, and the anisotropic val- 
ues are expected to be temperature independent for 
our S = 1/2 system over the temperature range of 
our measurements. Thus one expects the difference 
Xa{T) — xpiT) (a, P = a, h, c) to be nearly independent 
of temperature compared with the magnitude of either, 
where a subscript refers to the crystallographic axis along 
which the magnetic field is applied. 

To be more quantitative, we define the anisotropy in 
the intrinsic susceptibility as 



^Xaf}{T) = [Xa 



xS^](r) 



[X[3 



(73) 



which eliminates extrinsic anisotropy in the Curie- Weiss 
impurity contribution from the values calculated from 
the experimental data. For example, according to this 
definition, iS.Xac{T) is the difference between the upper- 
most and lowermost data sets in Fig. 39(b). The three 
Axa(3{T) anisotropies are plotted in Fig. It seems 
to us that the only reasonable explanation for the strong 
temperature-dependent anisotropies in Fig. ^ for two of 
the three data sets is that one or more of the x^ sus- 
ceptibilities is strongly temperature dependent near Tc, 
contrary to our initial expectations. Such a temperature 
dependence may be associated with the crystallographic 
and charge-ordering transitions which occur at or near 
the same temperature as the spin dimerization transi- 
tion, as discussed in the Introduction. 

One can make rather strong general statements about 
the magnetic susceptibility anisotropies and their tem- 
perature dependences as follows. Defining the Van Vleck 
susceptibility anisotropy Ax^J — x^ ~ X and simi- 
larly the spin susceptibility anisotropy ^Xiap^ = X^a™ ~ 
X^^™, from Eqs. ( |7l| ) one obtains an expression for the 
anisotropy l^XapiT) for a spin system in which the only 
anisotropy in ^'^p™ arises from anisotropy in the g factor, 
given by 



AxoMT) = AxXJ 



Ny^ligl 



J 



if) ■ 



The reduced spin susceptibility x*(^) is necessarily pos- 
itive, and it is isotropic for a Heisenberg spin system as 
noted above. Thus, if x^ and x^^ and therefore Ax^J 
are independent of temperature, the slope dl^Xap{T)/dT 
must have the same sign as the difference — 9^- 
As discussed in the next subsection, for NaV205, 
this difference has been reported to be positive for 
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FIG. 40. Temperature T dependences of the intrinsic 
anisotropy differences Axq/3 (q/3 = ac, be, ah) in the mag- 
netic susceptibilities along the a, b, and c axes in NaV205, 
as defined in Eq. (^3|). These data were obtained from the 
respective differences between the three pairs of data sets in 
Fig. 11(b). 

aP = ac and he and near zero for a(3 — ab, consistent 
with the slopes in Fig. However, in a simple ionic 
crystalline electric field model and with a positive spin- 
orbit coupling parameter for V one would predict that 
a Xq^ should increase with the negative deviation of ga 
from the free electron value g = 2. Thus, a particularly 
visible and puzzling discrepancy is that since (2 — ga) ~ 
(2 — gb) < (2 — gc) according to the reported ga values 

< 



below, on this basis one strongly expects x^^ 
X^^; thus two of the three Xq^ values should be about 
the same and smaller than the third one. Qualitatively 
contrary to this expectation, for T <^ Tc we observe in 

Eq. (g)thatxr«xr>xr- 

We will not emphasize or further discuss these puzzling 
discrepancies with expectation with respect to their pos- 
sible influence on ou r the oretica l modeling of our x{T) 
data in Sees. VIII B and VIII C, since at present there 
is no way to model, e.g., a temperature dependent Van 
Vleck susceptibility which changes rapidly near Tc, but 
the anisotropic susceptibility results and the above dis- 
cussion should be kept in mind. In the following two 
subsections the reported anisotropies in the g factor as 
measured using electron spin resonance (ESR) and in the 
Van Vleck susceptibility as deduced from nuclear mag- 
netic resonance (NMR) measurements will be discussed, 
respectively, in light of our anisotropic xC^) data. 

Anisotropy in the g faetor from ESR. Many 
ESR mciasurepafints have recently been reported for 
NaV205E3^n^El Each study found a signal with g ^ 2 
which was attributed to bulk V species, and the g values 
found in the, various studies were the same within the 
errors, e.g.,L3 



ffa- 56 = 1.972(2), 5^ = 1.938(2) 



(75) 



The powder-average value is g = \/ (ffa + 56 + 5c)/3 = 
1.961(2). The g values were found to be independent 
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of T down to 20 K, which is below Tc- From all these 
measurements, there is no indication that the S = 1/2 
Heisenberg Hamiltonian is not appropriate to the spin 
system in NaV205. Unfortunately, given the sensitivity 
of the ESR technique, we cannot be certain that these 
ESR results are representative of the bulk spin species 
in NaV205, because no quantitative measurements of 
the concentration of spin species observed in these mea- 
surements were reported. Although the (uncalibrated) 
ESR intensity versus temperature measurements approx- 
imately mirror the bulk susceptibility behavior in most 
(but not all) of these studies, it is still possible that the 
signal arises from a minority spin species that is cou- 
pled to the bulk spin system. An interesting related 
issue which has not been discussed in the literature is 
why the presumed bulk 5=1/2 species in NaV205 are 
observable to low temperatures T < 0.03 J/ A:b by ESR, 
where the AF exchange constant is ,///cb ^ 580 K (see 
below), whereas the bulk Cu+^ spins 1/2 in the high-Tc 
cuprates are not observable by ESR up to HOOK, which 
is « 0.2J/A:b where J/fce ~ 1600 K is only a factor of 2.8 
larger .□ 

In Ref. ^ the authors estimated the values us- 
ing the reported anisotropic g values obtained from 
ESR measurements, obtaining = Xb^ — 2.4 x 



10 °cm^/mol and % 



vv 



6.6 X 10 °cm^/mol, which 



were stated to be in agreement with the values from their 
K-x analysis discussed in the following subsection. These 
values do not agree with our T = values in Eq|_(72). 
In addition, from the values of Ohama et a/.J£il one 
obtains AxYa = ^X^b^ = 4.2 x 10^^ cm'^/mol, which are 
similar in magnitude but opposite in sign to our data in 
Eq. ([T^). If the strong change in each of Axac and Axfcc 
below Tc in Fig. ^ is due to a respective Ax^J which 
is strongly temperature dependent in this temperature 
range, an effect similar to that reported to occur from 
NMR measurements discussed in the next subsection, it 
is hard to understand why this change is not reflected in 
a distinct change in the reported temperature dependent 
anisotropy of the g values at if these g- value measure- 
ments are recording the characteristics of the bulk phase. 

Anisotropy in the Van Vleck susceptibility from NMR. 
From a so-called K-x analysis using NMR paramagnetic 
nuclear resonance shift K{T) data, combined with x{T) 
measurements, under certain assumptions x^^ can be 
obtained if K is proportional to X; with T as an im- 
plicit parameter. In this way, x^^ values have been 
obtained by Ohama and coworkers for NaV205 using 
23Na (Ref. |8l]) and ^^V (Ref. |l|) NMR measurements 
on the same aligned powder sample. The former ^^Na 
study yielded x^^ = 23 x 10~^cm^/mol below Tc and 
16 X 10~^ cm'^/mol above Tc, corresponding to a decrease 
of 7 X lO^^cm'^/mol at Tc. Their low temperature value 
is quite similar to our value in Eq. (|7^). 

The ^^V NMR study,Eil carried out above Tc, yielded 
2(1) X 10~^ cm'^/mol, roughly an order of mag- 
nitude smaller than obtained in the authors' first study 
(no comment was made about this discrepancy), and in 
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addition gave xT^ = 1(1) ^ 10~^cm'^/mol and x^^ = 
4(1) X 10~^cm^/mol. These values are significantly 
smaller than our values. We note that a K-x anal- 
ysis on the d} V"''^-, compound VO2 yielded x^^ = 
6.5 X lO-^cm^/molEa 



B. Modeling the susceptibility of NaV2 05 above Tc 



Turning now to the experimental x{T) data in Fig. pq , 
we have T"^^^ w 370 K. A ssuming the validity of the 
Hamiltonian (|l|), Eq. (30a) for the uniform Heisenberg 
chain yields the exchange constant J/k-Q 580 K. Then 
the gh value in Eq. ( [75| ) and our xo values at T = 
in Table ^ below, together with Eqs. (|l]) and (|7^), 
predict that the measured x"^^^ ^ 40 x lO^'^ cm'^/mol, 
which is similar to the measured values of « 44 and 
48 X 10~^cm'^/mol for the two crystals in Fig. re- 
spectively. We therefore proceeded to try to fit the data 
by the uniform chain model. The "Fif'^shown as the 
sohd curve in Fig. is a plot of Eqs. (|7l|), with x*(*) 
being the susceptibility of the uniform chain (Fit 2 above) 
and with the parameters 



Xo = 8 X 10-5 -, Cimp = 0, 

mol 



= 1.972, -^ = 580K. ("Fit") 
Kb 



(76) 



This "Fit" is not really a fit, since we just set the g 
and J values to those estimated above and then set xo 
so that the calculated curve is in the vicinity of the 
data, because no small change in the parameters can 
bring the theory in agreement with the data. It is clear 
that adjusting xo further will not improve the agree- 
ment, nor will including a nonzero impurity Curie con- 
stant Cimp. However, the shapes of the curve and the 
data are similar, so the agreement can be improved con- 
siderably (not shown) by simultaneously decreasing xo 
to « —10 X 10^5 cm'^/mol, which is not possible physi- 
cally according to Eqs. ( [tiI ) because it would require the 
Van Vleck susceptibility to be negative, and increasing 
g to the unphysically large value of « 2.4, while keep- 
ing J constant. These results are in disagreement with 
the conclusion of Isobe apd Ueda who found that the 
Bonner-Fisher predictionlU fitted their powder suscepti-. 
bility data from 50 to 700 K very well assuming g — 2£3 
We can only note that their x{T) data-iiai* not been 
quantitatively reproduced in either theirEllo or others' 
subsequent measurements on NaV2 05, including ours, 
and that the Bonner-Fisher prediction is not accurate at 
temperatures below ~ J/(4fcB) ~ 145 K as discussed in 
the Introduction. _ _ 

Lohmann et alVJi and Hemberger et aLLj also previ- 
ously concluded that the x{T) of NaV2 05 is not de- 
scribed (below 250 K) by the prediction for the S" = 1/2 
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Heisenberg chain- based on their fits by the Bonner- 
Fisher predictionEI to their x(T) deduced from ESR mea- 
surements up to 650 K. They suggested that additional 
exchange coupHngs may be required to explain the ob- 
served x{T)- We consider this possibility here by mod- 
eling the influence of possible interchain spin coupling. 
Because there are no accurate and generally applicable 
numerical calculations for this case that we are aware 
of, we utilize the following simple molecularpfield theory 
(MFT) prediction for the spin susceptibilityQL3 



1 



1 



z'J' 

~1~ 



(77) 



where XchairM-) reduced spin susceptibility of the 

isolated quantum S* = 1/2 uniform Heisenberg chain (our 
Fit 2 above). The parameter z' is the effective number of 
spins on other chains to which a spin in a given chain is 
coupled with effective (or average) exchange constant J'. 
To be consistent with our sign convention for the intra- 
chain exchange constant J, J' is positive for AF interac- 
tions and negative for ferromagnetic (FM) interactions. 
Equation ( [TT] ) is very accurate when \z' J / J\ <C iBtll 

We fitted the x{T) data above 50 K for the two crys- 
tals in Fig. |3^ by Eqs. ( pl[ ) and (|77|), where we fixed 
gb = 1.972 and Cimp = and allowed xo: J a-nd the prod- 
uct z' J' to vary. Very good fits were obtained, for which 
the fitting parameters are given in Table The fits are 
plotted as the solid curves in Fig. For the parame- 
ters of the two crystals taken together, the fitted J/fce 
= 584(9) K is the same as deduced above (580 K) from 
the temperature of the maximum in x{^^)^ a-nd the fitted 
Xo = 1.4(16) X 10"^cm'^/mol is similar to our results at 
low temperatures in Fig. ^(b). A moderately large and 
negative (FM) interchain coupling z'J'/J = —1.26(5) 
was obtained. This coupling is sufhciently strong that 
long-range magnetic ordering might be expected, but 
which is not observed, possibly due to magnetic frustra- 
tion effects. If the present mean-field interchain coupling 
analysis is correct, this interchain coupling should be ev- 
ident in the magnon dispersion relations observable by 
inelastic magnetic neutron scattering measurements. In- 
deed, moderately strong dispersions of 1.4 meV in each of 
two bands perpendicular to ±lie chains have in fact been 
observed by Yosihama et aZH in such measurements on 
single crystals. It remains to be seen whether the mag- 
nitude and sign of the interchain exchange coupling that 
we infer in the mean-field analysis are consistent with the 
dispersion relations deduced from the neutron scattering 
data. 

An alternative and/or additional mechanism which 
can produce a strong deviation of the measured x{T) 
of a uniform chain compound from that predicted for 
Heisenberg unifornx_and pajternating chains is the spin- 
phonon interaction.Eii£3'ETE£l At low T this interaction 
can lead to a spin-Peierls transition and can strongly 
modify x(Ti)-aJ3pve from that expected for the Heisen- 
berg chain.OEj Sandvik, Singh, and Campbell carried 
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FIG. 41. Fits of the magnetic susceptibility x vs temper- 
ature T from 50 K to 750 K for the two crystals in Fig. Kq 



with Cii 



0, 



9 



1.972 and the spin suscep- 



given by the molecular field theory (MFT) pre- 
fer coupled quantum 5 = 1/2 uniform Heisenber 



by Eqs. ([71| 
tibility 
diction (17 

chains. The fits are shown as the solid curves and the fit pa- 
rameters are given in Table The fits overlap the data so 
they are difficult to see; consequently they have been extrap- 
olated to higher and lower temperatures to show where they 
are. 

out a detailed QMC investigation of a spin-Peierls model 
in which the spin 1/2 interactions were modified by the 
presence of dynamical|4fluantum mechanical) dispersion- 
less Einstein phonons.Efl For particular values of the spin- 
phonon coupling constant and phonon frequency, they 
found that the effective exchange constant Jeff decreases 
strongly with increasing T, and at T = is about 27.3% 
larger than the bare J. Perhaps surprisingly, they found 
however that if the bare g factor is reduced by « 7% 
and the bare J by « 18% in the x*(t) predicted for the 
Heisenberg model, this model was then in good agree- 
ment with their QMC simulations for temperatures above 
Tc- A recent important extensive study of many finite- 
temperature properties of the same model usiiig-QMC 
simulations was carried out by Kiihne and L6w.e3 They 
found that for not too low temperatures, the suscepti- 
bilities for various Einstein phonon frequencies and spin- 
phonon coupling constants can all be scaled onto a uni- 
versal curve, given by that for the uniform Heisenberg 
chain, using only a suitably defined effective exchange 
constant Jeff > J. Contrary to the result of Ref. [5^, they 
found that a rescaling of the g factor was not necessary. 

TABLE V. Fit parameters for the magnetic susceptibility 
of two NaV2 05 crystals according to Eq. ( |7l| ) with Cimp = 0, 
g — 1.972 and the spin susceptibility x'^'"" given by the molec- 
ular field theory expression (7^) for coupled quantum S = 1/2 
uniform Heisenberg chains. 





Xo 


J/fcB 


z'J'/J 


crystal 


(10-5£h4) 

\ mol / 


(K) 




E083EF 


+2.8(2) 


577(2) 


-1.28(3) 


E106E 


-0.1(1) 


592(1) 


-1.23(2) 
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Our experimental results for NaV205 are not consistent 
with either of these theoretical studies, because as dis- 
cussed below Eq. (^) above, to force-fit the Heisenberg 
chain x(T) prediction onto the data requires an unphysi- 
cally large negative value of xo , as well as an unphysically 
large increase in g. 

On the other hand, our observed x{T) does not agree 
with the Heisenberg chain model (with a temperature- 
independent J), and in the next section we simul- 
taneously model the data both above and below Tc 
within the context of the Heisenberg chain model with a 
temperature-dependent J, where we find that J{T) above 
Tc is very similar in form to that deduced in the calcula- 
tions of Refs. and [5^. Thus it may be the case that the 
spin-phonon interaction is indeed important to determin- 
ing x{T) in NaV205, but where the effects on x{T) are 
somewhat different than calculated in the models. In par- 
ticular, the theoretical predictions may be substantially 
modified if phonon spectra appropriate to real materials 
were to be used in the calculations instead of dispersion- 
less Einstein phonons. 



C. Simultaneous modeling of the susceptibility of 
NaV2 05 below and above Tc 

Previous modeling of xi^) of NaV205 to extract the 
spin gap has usually been done at the lowest temper- 
atures without reference to the magnitude of x above 
Tc. Here we utilize our fit to the for the Heisen- 

berg chain to extract J above from the experimental 
data. Clearly, since the measured x{T) above Tc cannot 
be modeled within this framework using a temperature- 
independent J as shown in the previous section, it follows 
that if we are to remain within this framework, J, which 
is then evidently an effective exchange constant incorpo- 
rating additional physics of the material, must be tem- 
perature dependent. Then with J(T) fixed, we derive 
the T-dependent spin gap A(T) and exchange alterna- 
tion parameter S{T) near and below Tc directly from the 
measured x(T) data, which has not, to our knowledge, 
been carried out before for any system showing a spin- 
dimerization transition, using our x*(q!, t) fit function for 
the alternating chain. 

The specific procedure we adopted for modeling the 
Xb{T) measurement on each crystal consists of the fol- 
lowing six steps, where we fixed gb = 1.972 in steps 3j5. 



(1) The xiT) from 2 to 10 K is fitted by Eqs. (|7l|), 
setting ;^''P'" = because of the presence of the spin gap, 
thereby obtaining the parameters xoi Cimp, and 9. 

(2) Using these xo, Cimp and parameters, we solve for 
J(T) for T > 60K, or for T = 50 K only, using our "Fit 1" 
function for x*(i) of the Heisenberg chain, which is one 
end-point function of our x* {5, t) fit function, and fit the 
resulting J(T) by a polynomial in T for extrapolation 
below Tc; we used the extrapolation function J(T) = 
J{0)+aT^ + bT^. 
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FIG. 42. (a) Magnetic susceptibility x versus temperature 
T for three crystals of NaV205 in the low-T regime near 
the dimerization transition temperature Tc ~ 33-34 K. The 
crystal symbol designations are E083EF (•), E097A (o), ALl 
(filled squares). The solid curves are fits to the data below 



20 K by Eq. ( 71a ) , where the spin gap is assumed independent 
of T, and have been extrapolated to higher temperatures, (b) 
Magnetic spin susceptibility x^'""(T), obtained from the data 
in (a) by subtracting [xo + Ci mp / (T — 9)] appropriate to each 
crystal according to Eq. (71a). 



(3) With this J(T), or using J(50K) only, we fit- 
ted x{T) from 2 to 20 K, now including for 
the alternating-exchange chain [i.e., using our alternat- 
ing chain x* (<5, t) fit function] assuming a T-independent 
S (and A), and obtain a new set of xqj Cimp and 9 pa- 
rameters [in addition to d{0)]. 

(4) Steps 2 and 3 are repeated until convergence is 
achieved, which takes in practice only one additional it- 
eration. Note that we implicitly assume that xo, Cmp 
and 6 are independent of T, i.e., that the transition(s) at 
Tc do not affect them. 

(5) The experimentally determined molar spin suscep- 
tibility x^'""(^) is now computed by i nser ting the final 
X07 Cimp and 9 fit parameters into Eq. (71a). Then using 
the fitted J(T) or J(50K), the S{T) is computed using 
our fit function x* {S, t) for the alternating-exchange chain 
by finding the root for S, at each data point temperature 
T, of 
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FIG. 43. Exchange constant J versus temperature T for 
the three crystals of NaV205 in Fig. ^ in the temperature 
regime T > 60 K which is above the dimerization temperature 
Tc ~ 34 K. The solid curves are polynomial fits to the data 
between 60 and 150 K for the respective samples, which are 
extrapolated to T = as shown. 
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(6) In a separate step not associated with the fitting pro- 
cedure in steps 1-5, the spin gap A(r) is computed from 
6{T) determined in step 5 using an independently known 
function A^(5) = A{6)/J and our J(T) or J(50K). We 
used our A*(^) fit function in Eqs. (|63| ) for this purpose. 

We measured Xb{T) for nine different crystals from four 
different batches of NaV205 and now present illustrative 
results obtained in each of the above modeling steps 
2 to 4 (final iteration), 5 and 6 for three representative 
crystals. We will follow in graphical form the data mod- 
eling through successive steps for these three crystals to 
show how differences in one property between the crys- 
tals may or may not propagate through the next step(s) 



of the analysis, but we present the fitting parameters for 
all of the crystals in Table |v^ . 

The measured x{T) data below 50 K for the three crys- 
tals are shown in Fig. ^(a), where the fits below 20 K in 
step 4 are shown as the solid curves with parameters in 
Table VI. Crystals E097A and ALl are seen to have much 
lower levels of paramagnetic impurities than E083EF, as 
refiected in the impurity Curie constant, i.e., the magni- 
tude of the impurity Curie- Weiss upturn at low T. By 
subtracting the xo and impurity Curie- Weiss terms from 
the data, the spin susceptibility x^^"^{T) is obtained for 
each crystal as shown in Fig. ^(b). These data show 
good consistency below Tc for the three crystals, despite 
the differences in the xo values, the magnitudes of the 
Curie- Weiss impurity term and in the x{T) above Tc- 
The J(T) determined for the three crystals in step 2 are 
shown up to 300 K in Fig. J is found to decrease 
by ~ 10-20% upon increasing T from 60 to 300 K, which 
when T is scaled by J is sirailar to the fractional decrease 
predicted by Sandvik et alL3 due to the spin-phonon in- 
teraction. It is noteworthy that crystal E083EF, with by 
far the highest level of paramagnetic defects and/or im- 
purities, also has the largest J{T) and the largest change 
in J with T. 

Figures ^(a) and ^(b) show the spin dimerization pa- 
rameter S(T) and spin gap A(T) determined for each of 
the three crystals in the final modeling steps 5 and 6, 
respectively. Several features of these data are of note. 
First, there is a rather large variation in the dimeriza- 
tion parameter, d{0) = 0.028-0.040, between the three 
crystals, despite the fact that Tc = 33-34 K is nearly the 
same for the different crystals; the most impure crystal 
E083EF has the smallest i5(0), as might have been ex- 
pected. Despite this variability, these (5(0) values are all 
significantly smaller than the three values reported for 
various samples by different groups as determi ned using 
different techniques, which are listed in Table VII along 



TABLE VI. Fitted parameters in Eqs. (fn]), using the x*('^j't) fit function ( |5c| ) for the alternating-exchange chain, obtained 
by fitting the xiT) data in the range 2-20 K for nine crystals of NaV205, for an assumed g factor of 1.972. If an arror bar is 
not given for J(0), this value is J(50K) which was determined from a single data point near 50 K. The spin gap A(0) is not a 
fitted parameter, but is rather computed from the fitted alternation parameter 5(0) using Eqs. (^^. Similarly, the alternation 
parameter a(0) is computed from 5(0) using Eq. (141:). Note that all three measurements for crystal ALl were carried out in 
a field of 2 T, whereas all the other crystals were measured in a field of 1 T. 



Crystal 


XO 


C'imp 


-e 


J(0)/fcB 


m 


«(0) 


A(0)/fcB 




(10-^7^) 


(10-^ ^) 


(K) 


(K) 






(K) 


E082E 


6.82(4) 


1.123(5) 


0.46(1) 


710(4) 


0.0287(2) 


0.9442(4) 


101.1(10) 


E083B 


4.56(9) 


0.81(1) 


0.30(5) 


654(3) 


0.0327(3) 


0.9366(6) 


102.6(13) 


E083EF 


11.2(1) 


1.11(1) 


0.43(4) 


723(2) 


0.0279(4) 


0.9458(8) 


100.8(14) 


E083G 


6.55(3) 


1.112(3) 


0.45(1) 


688 


0.0298(1) 


0.9421(2) 


100.7(2) 


E083H 


4.24(3) 


0.780(3) 


0.33(1) 


650 


0.0332(2) 


0.9357(4) 


102.9(5) 


E083I 


4.18(6) 


0.946(7) 


0.32(2) 


657 


0.0329(4) 


0.9363(8) 


103.3(9) 


E097A 


5.92(3) 


0.170(3) 


0.25(4) 


589(2) 


0.0389(2) 


0.9251(4) 


104.6(6) 


E106E 


5.67(8) 


0.134(8) 


0.31(1) 


662(2) 


0.0332(5) 


0.9358(9) 


104.8(13) 


ALl iH\\a) 


12.22(6) 


0.221(7) 


0.46(6) 


598(1) 


0.0366(3) 


0.9294(6) 


101.5(8) 


ALl (iJ| 6) 


10.94(6) 


0.240(7) 


0.46(6) 


607(1) 


0.0352(3) 


0.9320(6) 


100.4(8) 


ALl Ih\\c) 


5.49(7) 


0.298(8) 


0.75(7) 


635(1) 


0.0337(4) 


0.9348(7) 


101.7(9) 
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FIG. 44. Alternation parameter 5 (a) and spin gap A (b) 
versus temperature T below 50 K for the three crystals of 
NaV205 in Fig. The nonzero S and A above the tran- 
sition temperature ~ 33-34 K are presumed to arise from 
spin dimerization fluctuations and concurrent spin gap fluc- 
tuations, respectively. 



lated information. 
On the other hand, the 



corresponding range of A(0)/fcB = 103(2) K for the three 
crystals is fractionally much smaller than that of (5(0). 
We infer that so me o f the discrepancies between the A(0) 
values in Table VIl reported for NaV205 by different 
groups may arise from differences in, e.g., the types of 
measurements which are used to determine A(0) and in 
the different analyses of those data, rather than from 
different A(0) values in the samples. The variability in 
i5(0) between the crystals in Fig. ^(a), compared with 
the lack of much variability in A(0) in Fig. 4j(b), evi- 
dently arises because S must be combined with J to ob- 
tain A, and the variations in the first two parameters 
must largely cancel. Thus, not surprisingly, the low-T 
^spm^y-j is governed by the spin gap A and not by S or 
J separately. 

The 6{T) data for our best crystals show very sharp, 
nearly vertical increases with decreasing T at Tc. We can- 
not extract a precise critical exponent /3 from our S{T) 
data due to the large temperature-dependent background 
above Tc, to be discussed shortly. However, rough fits be- 
low Tc by the expression S{T) ^ [l—T/TcY gave fi values 
consistent with the values (3 = 0.25(10) (Ref. |90| ) from in- 
frared reflectivity measurements, 0.34(8) (Ref. |9l|) from 
sound velocity measurements along the chain axis and 
0.35(8) (Ref. Efl) from thermal expansion measurements 
along that axis. We note that these values are a factor 
of two larger than the value of ~ 0.15 (Ref. |9^) inferred 
from x-ray diffuse scattering measurements. 

The data in Figs. |4^(a) and ^(b) clearly show the exis- 
tence of spin dimerization fluctuations and a spin pseudo- 
gap above Tc for each crystal, respectively, irrespective 
of the crystal quality as judging from the Curie- Weiss 
impurity term in the low-T x(^)i with magnitudes just 
above Tc of about 20% and 40% of 5(0) and A(0), 
respectively. This is a robust result, which was obtained 



TABLE VII. Exchange constant J, spin gap A, and alternation parameter 5 for NaV205 at the temperature(s) T as 
determined by the listed method for the sample with transition temperature Tc- The literature reference is given in the last 
column. Method abbreviations: magnetic susceptibility; Neutrons, neutron scattering; NMR, nuclear magnetic resonance; 
ESR, electron spin resonance; Cp, specific heat; Raman, Raman light scattering. 

Ref. 



sample 


Tc(K) 


J/fe (K) 


r(K) 


A/fcB (K) 


5 


Method 


Powder 


none 


529 


350 






X 


Powder 


33.9 


560 


35-700 






X 


Powder 


35.3 




7 


114 




Neutrons 


Aligned powder 






10-20 


98 




23Na NMR 


Crystals 


33 


441 


2-30 


85(15) 




X 


Crystal 


35 


560 


2-34 


92(20) 


0.10(2) 


ESR, X 


Crystal 


33.5 


560 


15-30 


100(2) 


0.107 


ESR 


Crystal 


34 


578 


250-650 


100 




ESR, X 


Crystal 


32.7,33.0 




1.8-12 


84(10) 






Crystal 






5 


88(2) 




Raman 


Crystal 


35 




10-35 


85(20) 




ESR 


Crystal 


34 


455 


15 




0.047 


Raman 


Crystal 






4.2 


94 




ESR 


Crystals 


33 




7-15 


67(5) 




Cp 


Aligned powder 


34.0 




11-20 


108 




''^V NMR 


Powder 


34.0 


491 


4-30 


77 




X 
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for each of the nine crystals we measured, which does 
not depend on the precise value of J [and resultant 
X''P'"(T, (5 = 0)] or the details of how J is determined 
above Tc, or even on the detailed formulation of the 
X*{oL^t) fit function for the alternating-exchange chain. 
For example, setting J to be a constant, equal to the 
value at 50 K, yields nearly the same A(r) near Tc as de- 
termined using a T-dependent J. Similarly, deleting the 
impurity Curie- Weiss term in the fit to the data above Tc 
changes the derived xo and J{T) or J(50K) somewhat 
as well as the detailed temperature dependence of the 
pseudogap A(r) above Tc but has little influence on the 
magnitude of A near Tc. Further, in a previous version of 
the QMC and TMRG x*(a, t) fit function (not otherwise 
discussed in this paper ), we did not enforce the require- 
ment (iii) in Sec. IV E that the transformed x*{^^T) sat- 
isfy dx* (<5, t) /d5\s=o = 0, and the same fluctuation effects 
above Tc were found using that fit function as using the 
present one, although these fluctuations were somewhat 
reduced in magnitude compared to the present results. 
Finally, these fluctuations are observable directly in the 



measured x(T) data in Fig. 42(a) as a rounding of the 



susceptibility curves above Tc. 

From Fig. the fluctuation effects persist up to high 
temperatures T > 50 K, although the fluctuation am- 
plitudes decrease with increasing T. Precursor effects 
above Tc have been reported in x-ray diffuse scatter- 
ing measiirementsLJ up to ^ 90 Ky, in. , uJ trasonic sound 
velocityEll and optical absorptionCJ'EScI measupements 
up to ~ 70 K, and in specific heat measurementsOO up 
to ^ 40-50 K, so it is not surprising that spin dimer- 
ization parameter fluctuations in Fig. ^(a), and a spin 
pseudogap in Fig. |4^(b) reflecting fluctuations in the spin 
gap, are found above Tc. 



D. Specific heat of NaV2 05 

In order to correlate the magnetic effects discussed 
above in NaV205 with thermal effects, we have carried 
out specific heat vs temperature Cp(T) measurements on 
the same crystal ALl, and a crystal E097 from the same 
batch as E097A, for which x{T) data were presented and 
modeled above. The results from 2K to 50 K for crystals 
E097 and ALl are shown in Fig. ^(a). Over this tem- 
perature range, the Cp(T) data for the two crystals agree 
extremely well, except in the range 33.0-34.2 K, i.e., in 
the vicinity of the transitions as will be discussed shortly. 
The shapes of the specific heat anomalies at Tc are not 
mean- field-like specific heat jumps as observed in, e.g., 
conventional superconductors, but instead are A-shaped 
anomalies. Thus, any attempt to define a (mean-field) 
"specific heat jump at Tc" is fraught with ambiguity. 
These shapes are retained in plots of Cp{T)/T vs T as 
shown in Fig. ^(a). This A shape has been observed 
previously, and variously attributed to fluctuation effects 
or a possible smeared-out first order transition. In view 
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FIG. 45. (a) Specific heat Cp vs temperature T for NaV205 
crystals E097 and ALl. (b) Expanded plots of the data in the 
vicinity of the transition temperatures of the two crystals, (c) 
Temperature derivative of xT vs T for the same crystal ALl 
as in (a) and (b) plus data for crystals E097A (from the same 
batch as E097) and E106E. The lines connecting the data 
points are guides to the eye. 

of the coupled structural, charge-ordering and spin 
dimerization transitions at Tc in NaV205 as discussed 
in the Introduction, their relative contributions to the 
specific heat anomalies are not clear, if indeed their con- 
tributions can be uniquely distinguished. 

Expanded plots of Cp{T) and Cp{T)/T versus T, 
shown in Figs, ^(b) and p^(b), respectively, reveal a 
sharp high peak at 33.4 K for crystal ALl, which is 
slightly split by O.IK in spite of the fact that the 
overall height of the anomaly is much larger than pre- 
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FIG. 46. (a) Specific heat Cp divided by temperature T 
versus T for NaV205 crystals E097 and ALL (b) Expanded 
plots of the data in the vicinity of the transition temperatures 
of the two crystals. The lines connecting the data points for 
crystal ALl are guides to the eye. 

viously reported for any crystal of NaV205. Two peaks 
are also observed for crystal E097, at 33.4 K and 33.8 K, 
which are more widely separated than for crystal ALl. 
From Fig. ^(b), the entropy under the anomaly (ies) for 
each crystal is about the same (see below). Comparing 
these results with the S{T) and A(T) data in Fig. ||, 
the larger splitting of the Cp{T) peak for crystal E097 
does not result in any major difference in the magnetic 
order parameter properties between the two crystals, al- 
though the transition onset is slightly rounded for crystal 
E097A compared to ALl. By using the Fisher relation,E2l 
d[x{T)T]/dT - C(T) where C{T) is the magnetic con- 
tribution to the specific heat, one obtains results which 
show the same features near Tc as does the specific heat, 
as shown in Fig. ^(c). Thus, careful scrutiny of the mag- 
netic properties can reveal the fine detail observed in the 
specific heat near Tc- In particular, this comparison sug- 
gests that both anomalies in the specific heat near Tc 
for each crystal are associated with and/or reflected by 
magnetic effects. 

The splitting of the transition into two apparent transi- 
tions that we report here was previously observed in ther- 
mal expansion, but not seen in their spoeific heat, mea- 
surements of a crystal by Koppen et a/.c2l The detailed 
origin of the transition splitting, and more fundamentally 



whether the splitting is instrinsic to ideal crystallograph- 
ically ordered NaV205, remain to be clarified. An essen- 
tial feature that any explanation must account for is that 
the temperature splitting between the two transitions in 
a crystal varies from crystal to crystal. 

Modeling. In this section we will only consider the 
model utilized above for analyzing our x(T) data, in 
which NaV205 consists, effectively, of isolated 5=1/2 
uniform or (below Tc) alternating-exchange Heisenberg 
chains, where the (average) exchange constant J shows, 
at most, only a smooth and relatively small change below 
Tc. For reasons which will become clear below, unfortu- 
nately we cannot use our specific heat data to extract 
detailed information about the magnetic subsystem in 
NaV205. However, other types of important information 
about the thermodynamics will be derived using various 
of the theoretical results presented and discussed previ- 
ously in this paper. i— in 

There have been two reportal3'E£l deriving the spin 
gap from Cp(T) data at T < 15 K. We first discuss 
the limits of this type of analysis. Using J{0)/kQ = 
600K,^) = 0.040, and A/kn = lOOK (see Table 0), 
Eqs. ( ^9b| ) and ( ^9c| ) predict that the magnetic specific 
heat C(T) in the dimerized phase at low temperatures 
T< (A//CB, Tc) is 



C(T) = 1.0 



100\3/2 



molKV T / 
1 ^ 3 / T y 
^ 100 ^ 4 VlOO/ 



-lOO/T 



(79) 



with T in units of K. Equation ( [79[ ) predicts that 
C(15K) = 0.026 J/molK, which is about 40 times 
smaller than the observed Cp(15K) « 1 J/molK (which 
must therefore be due to the lattice contribution) and 
hence is unresolvable at such low temperatures. Within 
this model, we must therefore conclude that the previ- 
ous estimates of the spin gap based upon modeling the 
low temperature specific heat were most likely artifacts 
of modeling the lattice specific heat. This can happen 
if one does not utilize the fact that the prefactor to the 
activated exponential term of the magnetic contribution 
C(T) is not an independently adjustable parameter, but 
is rather determined by the spin gap itself as we hav e 
previously demonstrated and emphasized in Sec. [I C 3| . 

A related question is whether the entropy associated 
with the transition(s) at Tc can be associated solely with 
the magnetic subsystem. The minimum possible esti- 
mate of the entropy of the transition is obtained from the 
Cp(T)/T vs T data in Fig. ||(b) by drawing a horizontal 
line from the data at the Cp(T)/T minimum at « 35. OK, 
just above Tc, to the data that the line intersects with 
below Tc at ~ 30.6 K, and then computing the area be- 
tween the line and the peak(s) above the line. In this 
way we obtain a value of 0.397 J/molK for crystal E097 
and 0.375 J/molK for crystal ALl. On the other hand, 
the maximum magnetic entropy of the S' = 1/2 uniform 
chain subsystem at Tc, using rough values J/k^ — 600 K 
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FIG. 47. Measured specific heat Cp vs temperature T up 
to 100 K for NaV205 crystal E097 (o). The solid curve is the 
background specific heat, which is the specific heat that would 
have been observed had no transitions or order parameter 
fluctuations occurred, determined as described in the text. 

and Tc = 34 K, is roughly S{Tc) « {2R/2,){k^Tj J) ^ 
0.31J/molK. Thus, the specific heat A anomaly at Tc 
cannot arise solely from the magnetic subsystem, since 
the minimum possible entropy of the transition is sig- 
nificantly greater than the maximum possible magnetic 
entropy at Tc- At the least, the remaining entropy must 
therefore be due to the crystallographic and/or charge- 
ordering transitions which occur at or close to the spin 
dimerization transition temperature as discussed in the 
Introduction. 

A potentially definitive and effective way to proceed 
from this point would be to quantitatively determine the 
magnetic contribution C(T) to the measured specific heat 
Cp(T) from x(T) at and near Tc, using a relationship 
between x(T) and C(T) such as the Fisher relation cited 
above, and then compare this result with Cp(T). From a 
comparison of Figs, ^(b) and 45(c), it seems clear that 
such a relation must exist, at least for temperatures near 
Tc, but the relationship between x{T) and C(T) near 
spin dimerization transitions has not yet been worked 
out theoretically. 

In the absence of such a formulation, we proceed to es- 
timate the change in the specific heat associated with the 
transition(s). In order to do this modeling, we must fit 
Cp(r) to higher temperatures than we have been dis- 
cussing so far. The Cp{T) data from 2 to 100 K for 
NaV2 05 crystal E097 are shown as the open circles in 
Fig. ^ As noted above, except in the immediate vicin- 
ity of Tc the Cp(r) data for crystal ALl are nearly iden- 
tical with those for crystal E097 up to at least 50 K, so it 
will suffice to model the data for crystal E097. The four 
modeling steps and the assumptions we employed are as 
follows. 

(1) We assume that critical and other order parame- 
ter fluctuations associated with the transition(s) at Tc 
make a negligible contribution to Cp(T) over some spec- 
ified high temperature (T ^ Tc) range. By subtract- 
ing the known magnetic contribution C(T) due to iso- 



lated Heisenberg chains [obtained using our fit function 
for C{kBT/J)] in this temperature range, we obtain the 
background lattice contribution C''^* (T) in the high tem- 
perature region. Also, since we have shown that C(T) is 
neghgible for T < 15 K, the measured Cp(T) in this T 
range is assumed to be identical to C''**(T) at these tem- 
peratures (we again neglect the possible but unknown 
specific heats associated with possible order parameter 
fluctuations in this range). Thus we obtain background 
lattice specific heats C'*^* (T) in high and low temperature 
ranges which are assumed unaffected by the transition(s) 
and associated order parameter fiuctuations. 

(2) We interpolate between the C^^^{T) determined in 
step 1 in the low- and high-temperature ranges to obtain, 
in the intermediate temperature range, what C''^'(T) 
would have been in the absence of the transition(s) and 
associated order parameter fluctuations. 

(3) We add the C(T) for isolated chains, used in step 1, 
back to the C^^^{T) derived in step 2 over the full tem- 
perature range of the measurements. This is the total 
background speciflc heat that would have occurred in the 
absence of the transition(s) and associated order param- 
eter fluctuations. Then we subtract the total calculated 
background speciflc heat from the measured Cp(T) data. 
This difference AC(T) should hopefully be a reasonable 
estimate of the change in the speciflc heat associated with 
the transition and order parameter fluctuations, includ- 
ing all lattice, charge and spin contributions. AC(T) 
must go to zero, by construction, at temperatures above 
the lower end of the high temperature region fltted in 
step 2. 

(4) Finally we integrate AC/T with T up to and be- 
yond Tc to obtain the change in entropy AS{T) asso- 
ciated with the transition and order parameter fluctua- 
tions. AiS'(T) must become constant, by construction, at 
temperatures above the lower end of the high tempera- 
ture region fitted in step 2. 

In the following we will present and discuss the results 
in each of the four steps of our modeling program de- 
scribed above. 

Step 1. Here we fitst use our C(T) fit function for the 
numerical C(T) data,El which was given in Eqs. (Q), to 
extract C'''*(T) in the high-temperature region above Tc. 
For consistency with our analysis of the susceptibility in 
Sec. VIII C , we use the temperature-dependent J(T) de- 
rived in that section for crystal E097A (see Fig. ^) when 
computing C(T). The background C(T) thus estimated 
for crystal E097, i.e., the values which would have been 
observed if no transition(s) at Tc or associated order pa- 
rameter fluctuations had occurred, is shown in Fig. |4q . 
Comparison of these data with the measured Cp(T) data 
in Figs. ^ and ^ shows that this C(T) is a small, but 
non-neghgible (^ 1%), fraction of Cp(T) above Tc. On 
the other hand, C(T) is much larger than the observed 
Cp(T) at low temperatures, because in this temperature 
range C cx T whereas Cp(T) = C^'^^T) oc T^. The 
C^^\T) = Cp{T) - C(T) in the high temperature 
(60-100 K) region is shown in Fig. E^, together with 
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FIG. 48. Magnetic specific heat C vs temperature T for 
NaV2 05 crystal E097, calculated for uniform chains with ex- 
change constant J{T) determined from the analysis of the 
susceptibility data for crystal E097A. The values are those 
which are assumed to have been observed had no transitions 
or associated order parameter fluctuations occurred. 
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FIG. 49. Background lattice specific heat C ^' vs temper- 
ature T for NaV205 crystal E097 (o). The data shown, in 
the temperature ranges 2-15 K and 60-101 K, were fitted by 
a polynomial; this interpolation fit is shown as the solid curve. 
This background, including the curve in the interpolated tem- 
perature region, is the lattice specific heat assumed to have 
been observed had no transitions or associated order param- 
eter fiuctuations occurred. 

C^^\T) = Cp(T) in the low temperature (2-15 K) region. 

Step 2. In this step we must interpolate C''**(T) be- 
tween the low- and high-temperature regions, i.e., in a 
broad temperature range spanning the transition region. 
The best way to do this would be to determine C'^'(T) 
directly from Cp(T) measurements on a suitably chosen 
reference compound, but such measurements have not 
yet been done. At first sight, a physically realistic possi- 
bility might be to interpolate the low and high tempera- 
ture C''^* (T) data using the Debye specific heat function; 
however, this method is questionable because the Debye 
temperature 8d in real materials can be rather strongly 



temperature dependent within the temperature range of 
interest here. The Debye function for the molar lattiea 
specific heat at constant volume C^'^^y'^ (T) is given byE3 



(^Dcbyc(y) = 9rR 



7^ \ 3 r^r>/T 



(od) Jo 



(e^ 



T) 



■ dx 



(80) 



where r is the number of atoms per formula unit (?' — 8 
here) and R is the molar gas constant. We attempted 
to fit our C'^'(T) data for the temperature ranges 2- 
15 K and 40-100 K to 80-100 K by Eq. (|^). The fits 
parametrized the data very poorly. We obtained a more 
reasonable fit by allowing r to be a fitting parameter, 
yielding a fitted value r ~ 4, but the data were still poorly 
fitted, due to too much curvature in the Debye function 
in the high temperature region. Therefore, we were led 
to interpolating between the low- and high-temperature 
regions using a simple polynomial interpolation function. 

To obtain the background lattice specific heat interpo- 
lation function, we fitted the combined C'''*(T) data (a 
total of 141 data points) in the low and high temperature 
ranges 2-15 K and 60-101 K, respectively, by polynomials 
of the form 



E 

n— 3 



Tn 



(81) 



The minimum summation index n = 3 is set by the ex- 
pected Debye low-temperature behavior of the lattice 
specific heat. The maximum value n^'^^ was varied to 
see how the fit parameters and variance changed. In ad- 
dition, for checking the final fits we fitted the C''''(T) 
data in the 2-15 K low-T range together with C'^^'^T) 
data in a high-temperature range varying from 40-101 K 
to 90-101 K. We found that the most stable fits were 
for 71™^'^ = 7 and 8. For both values, the fit did not 
visibly change when the lower limit of the upper tem- 
perature range of the fitted data was varied from 60 to 
70 K. We chose to use the fit for n"^'^^ = 7 because in 
this case the fit was also stable for lower limits of 50 and 
80 K. This stability allows one to be confident that the 
interpolation of the fit between the fitted low- and high- 
temperature ranges is an accurate representation of the 
background lattice specific heat in the interpolated inter- 
mediate temperature range. The fit for the temperature 
ranges 2-15 K and 60-101 K is shown as the solid curve 
in Fig. ^ The absolute rms deviation of this fit from 
the fitted data is quite small, (Tims = 0.046 J/molK. The 
curve over the full temperature range 2-101 K represents 
the background lattice specific heat C''*'(T) expected in 
the absence of any transitions or order parameter fluctu- 
ations. 

Step 3. Adding the magnetic background specific 
heat contribution C{T) obtained in step 1 to the lattice 
background specific heat contribution C^^'^{T) obtained 
in step 2 gives the total background specific heat, which 
is plotted as the solid curve in Fig. |4^. We reiterate that 
this background is interpreted as the specific heat that 
would have been observed had the transition(s) at Tc 
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FIG. 50. Temperature T dependence of the change in the 
specific heat AC (a) and in the entropy AS (b) in NaV2 05 
crystal E097 due to the transition(s) at Tc ~ 34 K as well as to 
crystallographic, magnetic and charge order parameter fluc- 
tuations associated with this (these) transition (s). The occur- 
rances of negative AC and AS values at low temperatures are 
real effects due to loss of magnetic specific heat and magnetic 
entropy, respectively, at these temperatures due to the open- 
ing of the spin gap at Tc. By construction, AC(T > 60 K) =0 
and AS(T > 60K) = constant. The actual order parameter 
fluctuation effects likely extend to temperatures higher than 
60 K. 

and associated order parameter fluctuations not oc- 
curred. The difference AC between the measured Cp(T) 
and the total background specific heat is plotted versus 
temperature in Fig. |o|(a). As would have been qualita- 
tively anticipated, AC is negative below about 16 K due 
to the loss of magnetic specific heat at low temperatures 
arising from the opening of the spin gap at Tc . This nega- 
tive AC does not arise from a problem in our polynomial 
interpolation C'''*(T) fit function or from our C(T) func- 
tion; these functions are both positive for all T > 0. Since 
the magnetic background contribution is proportional to 
T and the lattice background contribution [which is as- 
sumed not to change below 15 K due to the occurrence of 
the transition(s)] is proportional to at low T, open- 
ing a spin gap at Tc must necessarily lead to a negative 
AC at sufficiently low temperatures since the magnetic 
contribution then becomes exponentially small there. 



Step 4- Finally, we can compute the change AS in the 
total entropy of the system versus temperature due to the 
transition(s) and associated order parameter fluctuations 
by integrating AC(T) from step 3 according to A5(T) — 

J^[AC{T)/T]dT. The result is shown in Fig. |^(b). The 
entropy change is negative below about 22 K, due to the 
loss of magnetic entropy at low temperatures associated 
with the loss of magnetic specific heat as just discussed. 
From conservation of magnetic entropy, this lost entropy 
must reappear at higher temperatures. 

By construction, step 2 requires that AC(T > 60 K) 
= and consequently A5(T > 60 K) = const. This 
requirement is not desirable, but we had to enforce it 
to ensure that the AC(T) and A5'(T) derived at lower 
temperatures were accurate. Since the effects of the or- 
der parameter fluctuations are likely to continue to be 
present at temperatures higher than 60 K, the AC(T) 
and AS'(T) at temperatures at and near 60 K in Fig. |5^ 
are lower limits. 

The net change in the entropy at 60 K in Fig. ^(b) due 
to the occurrence of the transition(s) at Tc ~ 34 K and 
associated order parameter fluctuations above and below 
Tc is AS'(60K) = 2.28J/molK. This is far larger than 
the maximum possible change AS^^^ = 0.556 J/mol K 
in the magnetic entropy at this temperature obtained 
from Fig. where this value is just the maximum pos- 
sible entropy of the magnetic subsystem at this temper- 
ature, confirming our qualitative conclusion above based 
on very rough arguments. In particular, our quanti- 
tative analysis indicates that at least 76% of the en- 
tropy change at 60 K must arise from the lattice and 
charge degrees of freedom, and only a minor fraction 
(< 24%) from the magnetic degrees of freedom. Simi- 
larly, at Tc = 33.7K, we obtain AS* = 1.38 J/molK and 
AS^,^^ = 0.311J/molK, yielding AS^,^^/AS < 23% at 
Tc. 

As a closing remark for this section, it is clear from 
Fig. ^ and the discussion in the above two paragraphs 
that AC and AS do not saturate to their respective 
high temperature limiting values until a temperature of 
at least 60 K is reached, which is almost twice Tc. The 
present analysis of the thermal behavior of NaV2 05 thus 
lends strong support to our independent analysis and in- 
terpretation of our magnetic susceptibility data for this 



compound in Sec. VIII C 



IX. 



SUMMARY AND CONCLUDING 
DISCUSSION 



We have shown that the high-accuracy numerical 
Bethe ansatz calculations of the magnetic susceptibil- 
ity x*{t) for the S = %J2 uniform Heisenberg chain by 
Kliimper and Johnston 1 are in excellent agreement with 
the theory of LukyanovQ over 18 decades of temperature 
at low temperatures. An independent high precision em- 
pirical fit to these data was obtained over 25 decades of 
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temperature which we found useful to determine the ac- 
curacy of our TMRG x*(i) calculations. The magnetic 
specific heat dataEi for the uniform chain at very low tem- 
peratures was also compared with the theoretical predic- 
tions of Lukyanov, and extremely good agreement was 
found over many decades in temperature. We formulated 
an empirical fit function for these data which is highly 
accurate over a temperature range spanning 25 orders of 
magnitude; the infinite temperature entropy calculated 
using this fit function is within 8 parts in 10® of the ex- 
act value. We used both of the above fit functions to 
model our respective experimental data for NaV2 05 in 
later sections of the paper. We expect that they will be 
useful to other theorists and experimentalists as well. 

We have carried out extensive QMC simulations and 
TMRG calculations of x*(a,t) for the spin 5=1/2 
antiferromagnetic alternating-exchange Heisenberg chain 
for reduced temperatures t = k^T/ Ji from 0.002 to 10 
and alternation parameters a = Jij J\ from 0.05 to 1, 
where J\ (J2) is the larger (smaller) of the two alter- 
nating exchange constants. An accurate global two- 
dimensional (a, fit to these combined data was ob- 
tained, constrained by the fitting parameters for the ac- 
curately known x*(^) for the a parameter end points, the 
dimer (a — 0) and the uniform chain (a — 1), resulting in 
an accurate fit function over the entire range < a < 1 of 
the alternation parameter. Our fit function incorporates 
the first four terms of the exact high-temperature series 
expansion in powers of 1/f, which allows accurate extrap- 
olation to arbitrarily high temperatures. This function 
should prove useful for many applications including the 
modeling of experimental x(T) data as we have shown. 

Our x*(a,t) fit function for the alternating chain can 
be easily transformed (as we have done) into an equiv- 
alent fit function x*((5, i) in the two variables 8 = (Ji — 
J-i)l{^J) and t = ksT/J, where the average exchange 
constant is J = (Ji -I- J2)/2. This is a more appropri- 
ate function for analyzing experimental xC^) data for 
5 = 1/2 Heisenberg chain compounds showing dimeriza- 
tion transitions (such as a spin-Peierls transition) which 
result in an alternating-exchange chain with a small value 
of 6 at low temperatures. Once J has been determined 
by fitting our function for 6 = to the experimentally 
determined spin susceptibility x'*^™(^) data above the 
transition temperature, the alternation parameter S is 
uniquely determined by our fit function at each temper- 
ature below the transition temperature from the value of 
-^spin ^Yig^^ temperature. One can then find the spin 
gap A(r) using an independently known A* ((5). 

Our QMC and TMRG data and fit for x (a, t) are in 
good agreement with previous calculations based on ex- 
act diagonalization of the nearest neighbor Heisenberg 
Hamiltonian for short chains with a = 0.2, 0.4, 0.6, 
0.7, and 0.8, extrapolated to the thermodynamic limit, 
by Barnes and Riera.D However, the numerical and ana- 
lytical theoretical predictions of Bulaevskiijj which have 
been used extensively in the past by experimentalists to 
model their x(T) data for weakly-dimerized chain com- 



pounds, are found to be in poor agreement with our re- 
sults and should be abandoned for such use in favor of otnj 
fit function. Similarly, the previously used fit functionEa 
for the Bonner-Fisher calculation of x* (t) for the uniform 
chain {a — 1) should be replaced by one of our two. fit 
functions for the most accurate calculation to dateO of 
X*{t) for the uniform chain. 

An important theoretical issue in the study of the alter- 
nating exchange chain is how the spin gap A* (5) evolves 
as the uniform chain limit is approached (5 — *■ 0, a — *■ 1). 
We formulated a fit function for the temperature depen- 
dence of our TMRG susceptibility x*(«j t) calculations at 
low temperatures, which was used to extract the depen- 
dence A*{S) in this regime. We find that the asymptotic 
critical regime is not entered until, at least, S ^ 0.005 
(a }^ 0.99). We compared our spin gap data with many 
literature data. We formulated a fit function for our 
spin gapj-data together with those of Barnes, Riera, and 
TennantEj which quite accurately covers the entire range 
0<S <1. 

In the remainder of this paper, we showed how the 
above theoretical results could be used to obtain detailed 
information about real systems. As a specific illustration, 
we carried out a detailed analysis of our experimental 
x{T) and specific heat Cp(T) data for NaV205 crystals. 
This compound shows a transition to a spin dimerized 
state below the transition temperature Tc « 34 K. We 
used one of our two %*(<) fit functions for the uniform 
Heisenberg chain to model the xiT) above Tc, where we 
found that the experimental x(^) is not in quantitative 
agreement with the prediction for the uniform Heisenberg 
chain. A model incorporating a mean-field ferromagnetic 
interchain coupling between quantum 5=1/2 Heisen- 
berg chains fits the experimental data very well with rea- 
sonable parameters. It remains to be seen whether the 
inelastic neutron scattering measurements of the magnon 
dispersion relationsEJ are consistent with our derived in- 
trachain and interchain exchange constants. 

In an alternate description, we modeled the deviation 
in the measured x(^) of NaV205 above 60K > Tc from 
the Heisenberg chain model (with fixed exchange con- 
stant J) as due to a temperature-dependent J. We found 
that this J decreases with increasing T up to 300 K in 
a manner very similar.. to Jcg{T) predicted by_Sandvik, 
Singh and CampbelO and Kiihne and L6wE£| for the 
spin-Peierls chain. Our J(T) cannot however be com- 
pared directly with their Jeg{T) because the two quan- 
tities are defined differently. They found that by defin- 
ing an appropriate effective exchange constant Jeff, their 
resulting susceptibility xi^BT / Jcs) is universal at the 
higher temperatures for various Einstein phonon frequen- 
cies and spin-phonon coupling constants. This function 
agrees well with the xi^BT/J) for the 5 = 1/2 AF uni- 
form Heisenberg chain at these temperatures. As we 
discussed, these x(^) calculations are not applicable to 
NaV205, possibly because the calculations do not incor- 
porate realistic phonon spectra. 

Below Tc, we used the J(T) extrapolated from above 



50 



60 K and our global fit function for the alter- 

nating Heisenberg chain to determine the temperature- 
dependent alternation parameter S{T), and then the spin 
gap A(r) from S{T), directly from the x{T) data. We 
find that the A(0)/fcB values for nine single crystals of 
NaV205 are in the range 103(2) K. This result is in 
agreement, within the errors, with many previous anal- 
yses of data from various types of measurements for 
this compound by other groups. However, our values 
of 6{0) = 0.034(6) for various crystals are significantly 
smaller than previous estimates. We not e that the two 
estimates with ^(0) ~ OA in Table VI] were obtained 
using Bulaevskii's theoryQ for the alternating-exchange 
chain, which we have shown is not accurate at low tem- 
peratures in the relevant alternation parameter range. 

The dispersion of two one-magnon branches perpen- 
dicular to the chains observed in the neutron scattering 
measurements has been recently explained quantitatively 
by Gros and Valenti assumingj-that a zig-zag charge or- 
dering transition occurs at Tc.M They also predict that 
(5(0) ~ 0.034. This is within our range of (5(0) values in 
spite of the fact that we assumed that J{T) is either con- 
stant or increases slightly with decreasing T below Tc, 
contrary to their prediction that J decreases below T^. 
Gros and Valenti made no predictions for xi^), i5(T), 
A(T) or C{T), so comparisons with our results for these 
quantities are not pcisaible. We note that Kliimper, Rau- 
pach, and SchonfeldEj obtained a good fit to the xiT) 
data below Tc for the spin-Peierls compound CuGeOa 
within the context of a spin-Peierls model containing 
frustrating second-neighbor interactions and static spin- 
phonon coupling. 

We discovered that A(T) [and S{T)] of NaVaOs does 
not go to zero at Tc, indicating the existence of a spin 
pseudogap above Tc with a large magnitude just above 
Tc of « 40% of A(0); the pseudogap is present up to 
at least 50 K with a magnitude decreasing with increas- 
ing T above Tc- To our knowledge, this pseudogap has 
not been reported previously, and there are as yet no 
theoretical predictions for the magnitude or tempera- 
ture dependence of this pseudogap. The pseudogap is 
strongly reminiscent of the spin pseudogap derived by 
one of us using x(^) measurements above the transition 
temperature of inorganicj—nuasi-one-dimensional charge 
density wave compounda£3 as predicted theoretically by 
Lee, Rice, and AndersonO long before those observations 
were made. Similar to that case, in the present system 
one may think of the pseudogap as the rms fluctuation 
in the spin gap above Tc, with an associated reduction 
in the magnon density of states at low energy. In this 
interpretation, the pseudogap in NaVaOs should be ob- 
servable in high resolution quasielastic neutron scatter- 
ing and other spectroscopic measurements probing the 
low energy magnetic excitations. 

Finally, we carried out an extensive modeling study 
of our specific heat data for NaVaOs crystals, using the 
same model that we used to analyze our susceptibility 
data. The most important part of this study is that we 



have been able to determine a limit on the relative con- 
tributions of the magnetic and lattice/charge degrees of 
freedom to the entropy associated with the transition(s) 
at Tc- We find that at least 77% of the change in the 
entropy at Tc must arise from the lattice and/or charge 
degrees of freedom, to which the spin degrees of freedom 
must of course be coupled, and that the spin degrees of 
freedom themselves contribute less than 23 % of this en- 
tropy change. Our results also indicate that order param- 
eter fluctuation effects are important in the specific heat 
up to at least 60 K, strongly confirming the above simi- 
lar and independent conclusion based on our modeling of 
our magnetic susceptibility data for the same crystals. 
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